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ABSTRACT 

The  size  and  cost  of  microcomputers  continue  to  decrease 
while  their  memory  capacity  and  execution  speed  increase. 
These  advances  should  result  in  small,  inexpensive  machines 
attaining  the  same  computing  power  as  current  mainframe 
models.   The  interim  need  is  to  adapt  general  finite  element 
codes  to  present  day,  less  capable  microcomputers.   This 
thesis  explores  the  program  structure,   memory  management, 
I/O  procedures  and  equation  solving  methods  necessary  to 
accomplish  that  task.   The  equation  solving  capacity  and 
speed  of  the  Apple-II  Plus  Personal  Computer  Systems  and  the 
Hewlett-Packard  System  45(A)  Desktop  Computer  are  compared. 
A  finite  element  program  for  the  static  analysis  of  space 
trusses  is  presented,  as  adapted  to  and  tested  on  the  Apple-II 
Plus.   The  program  output  may  be  printed  in  either  English  or 
Fr ench . 
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I  .   INTRODUCTION 

A.   PURPOSE  AND  SCOPE  OF  THE  INVESTIGATION 

In  the  preface  of  the  third  edition  of  his  classic  refer- 
ence on  the  finite  element  method,  Zienkiewicz  [1]  made  an 
assessment  of  the  "state-of-the-art"  of  this  most  useful  analy^ 
tical  technique.   In  1977,  Zienkiewicz  described  a  discipline 
which  had  ammassed  a  collection  of  nearly  8000  published  ref- 
erence works  in  an  ever  widening  variety  of  applications, 
and  which  had  evolved  far  beyond  its  roots  as  an  engineering 
tool  for  linear  analysis  of  solid  mechanics  structures. 
Present  applications  commonly  include  heat  transfer,  fluid 
mechanics  and  non-linear   solid  mechanics,  as  well  as  the  com- 
plete static  and  dynamic  analysis  of  engineering  structures. 

A  more  geneal  statement  may  be  made  [1]  in  describing  the 
finite  element  method  as  "...  a  general  discretization  pro- 
cedure of  continuum  problems,  posed  by  mathematically  defined 
statements."   The  method  now  becomes  the  means  of  approximat- 
ing a  continuous  problem  with  a  discrete  model  composed  of 
finite  elements,  as  opposed  to  "infinitesimal"  elements  in 
the  sense  of  the  calculus.   In  some  cases,  in  the  limit,  as 
the  number  of  elements  increases  to  infinity  and  the  physical 
dimensions  these  elements  represent  decrease,  the  solution 
of  the  model  becomes  exact.   It  is  this  generality  of  method 
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and  application  which  has  been  one  of  the  two  most  important 
factors  influencing  the  widespread  use  of  the  method  in 
scientific  research  and  industrial  applications. 

The  second  factor  influencing  the  spread  of  the  finite 
element  method  has  been  the  availability  of  high  speed,  high 
capacity,  general  purpose  digital  computers.   It  is  well 
known  that  the  method  remained  in  its  infancy  until  general 
purpose  digital  computers  became  available  in  reasonable 
numbers.   Today  technological  advances  in  microelectronics 
and  solid  state  devices  have  drastically  reduced  the  real 
cost  of  computer  time.   Advances  in  the  art  of  building 
static,  random  access  memory  now  allow  the  use  of  Direct 
Memory  Access  techniques  in  small  (16+  bit  word)  desktop 
machines.   When  these  machines  are  interfaced  with  magnetic 
hard  disk  drives,  a  small  user  easily  has  the  potential  of 
16  Megabytes  (8  binary  bits  equals  one  byte)  of  high  speed 
storage,  at  a  small  fraction  of  the  cost  of  a  large  main- 
frame computer  [2].   These  developments  are  truely  revolu- 
tionary.  Very  "friendly"  and  general  structural  finite 
element  codes,  some  with  appealing  interactive  3raphics 
capabilities,  may  soon  become  available  to  the  smallest  of 
organizations . 

In  addition  to  the  new  machines  now  on  the  horizon,  there 
are  a  very  large  number  of  smaller,  less  capable  computers 
presently  in  service.   These  smaller  machines,  often  des- 
cribed as  "microcomputers,"  vary  in  word  length  and  storage 
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capacity,  generally  having  from  16  to  256  thousand  bytes  of 
random  access  memory.   They  support  a  wide  variety  of  periph- 
eral equipments,  the  most  indispens ib le  of  which  are  high 
speed,  magnetic  disk  drives.    These  smaller  machines  are  pro- 
portionately less  expensive  than  larger  computers,  which 
accounts  for  the  very  great  demand  for  these  systems  in  industry, 
government  and  the  private  sector.   Many  of  the  most  popular 
models  began  as  "personal"  or  leisure  computers  for  hobbiests, 
but  rapidly  moved  into  small  business,  medical  and  financial 
management  applications  as  the  utility  of  the  machines  was 
realized  and  software  became  available. 

If  one  looks,  however,  for  finite  element  software  or 
attempts  to  find  evidence  of  this  widespread  analytical  tech- 
nique being  implemented  in  any  substantial  manner  on  micro- 
computers, the  results  are  disappointing.   It  is  the  thesis 
of  this  investigation  that: 

(1)  Substantial  problems  in  finite  element  analysis  can 
be  solved  on  microcomputers. 

(2)  Microcomputers  can  solve  substantial  finite  element 
problems  with  acceptable  accuracy  and  to  reasonable  precision. 


Disk  storage  extends  the  usefulness  of  microcomputers  by 
avoiding  or  softening  the  impact  of  some  of  the  limits 
imposed  by  having  too  little  main  "core"  memory. 
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(3)  Although  the  central  processor  clock  speed  of  most 
popular  microcomputers  is  much  less  than  that  of  a  mainframe 
computer,  and  in  spite  of  the  interpretive  nature  of  most 
microcomputer  languages,  substantial  finite  element  problems 
can  be  solved  on  these  machines  in  an  acceptable  amount  of 

t  ime  . 

(4)  General  guidelines  can  be  established  for  using  inno- 
vative program  structuring,  memory  management,  software  overlay 
schemes,  data  base  management  and  high  speed,  disk  storage  to 
allow  the  implementation  of  substantial  finite  element  codes 
on  small  popular  microcomputers. 

B.   CHOICE  OF  MACHINES  AND  PERIPHERAL  DEVICES 

The  choice  of  machines  was  influenced  most  heavily  by  the 
types  of  systems  available  at  the  Naval  Postgraduate  School. 
A  Hewlett-Packard  System  45  Desktop  Computer  was  available 
within  the  Mechanical  Engineering  Department  and  was  a  logical 
first  choice.   Previous  work  had  already  been  done  on  imple- 
menting a  finite  element  code  on  this  system  [3]  and  the 
potential  existed  for  increasing  the  capabilities  of  the 
program  by  adding  additional  finite  elements,  changing  from 
an  "in-core"  to  an  "out -o f -cor e"  equation  solver,  and  using 
magnetic  floppy  disk  instead  of  magnetic  cartridge  tape  mass 
storage.   These  systems  are  not  uncommon  within  the  Department 
of  Defense  and  the  Navy  and  the  potential  benefit  from  an 
improved  program  was  judged  to  be  worth  the  additional  effort 
vhich  would  be  required  to  upgrade  the  code. 
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The  System  45  is,  however,  nearly  outside  the  boundaries 
of  the  group  of  machines  previously  described  as  "popular 
microcomputers."   The  base  price  of  the  computer  alone,  ex- 
clusive of  the  disk  drives,  was  over  $29,000.00  in  1977.   A 
new  System  45  (HP  9845C),  upgraded  and  having  slightly  more 
computing  capability,  would  have  cost  about  $39,500.00  if 
purchased  in  1980  [4]  and  still  be  without  disk  mass  storage. 
It  was  decided,  therefore,  to  also  choose  a  system  which 

would  represent  those  on  the  lower  end  of  the  cost  spectrum. 

2 
Several  Apple-II  Plus   Personal  Computer  Systems  were  avail- 
able, having  a  total  price,  including  multiple  disk  drives, 
of  approximately  $6,000.00.   In  addition  to  having  a  low  cost, 
the  Apple-II  Plus  systems  also  had  the  capability  of  the 
FORTRAN  Language  Option,  not  available  on  the  HP  9845A. 

For  the  reasons  enumerated  above,  it  was  decided  to  in- 
vestigate the  equation  solving  capability  of  the  Apple-II 
Plus  Personal  Computer  and  the  Hewlett-Packard  System  45 
Desktop  Computer  (HP  9845A)  and  to  attempt,  should  adequate 
time  be  available,  the  installation  of  representative  finite 
element  codes  on  both  systems. 

1  .   The  Hewlett  Packard  9845  Desktop  Computer 

The  System  45  was,  in  1980,  the  t op-of - the-1 ine  model 
of  Hewlett-Packard  desktop  computers.   As  with  all  of  this 


"Apple-II  Plus  is  a  registered  trademark  of  the  Apple 
Computer  Company. 
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firm's  computers,  technical  details  about  the  hardware  are 
closely  held  by  the  company  and  few  are  available  in  users 
manuals  supplied  with  the  machine  or  in  the  Hewlett-Packard 
Journal.   References  [5]  and  [6]  and  the  installation  docu- 
ments for  the  system  provided  the  information  necessary  to 
conduct  the  investigation. 

Only  a  very  brief  summary  of  the  most  important  details 
appears  in  the  following  sections,  and  readers  intending  to  use 
programs  presented  in  this  thesis  must  become  familiar  with 
the  machine  in  general  and  their  own  systems  in  particular. 

a.   System  Configuration 

The  system  configuration  during  the  investigation 
is  shown  in  the  table  below.   Only  the  minimum  devices  and 
options  which  are  required  to  run  the  program  presented  in 
this  thesis  are  shown.   Note  that  the  possession  of  a  more 
capable  system  than  the  one  shown  (e.g.  one  hard  disk  instead 
of  dual  flexible  disk  drives)  will  not  preclude  use  of  the 
proposed  software.   The  reader  need  only  change  the  "mass 
storage  unit  specifier"  [6]  and  "select"  [5]  codes  to  match 
his/her  particular  configuration. 

Although,  under  the  "unified  mass  storage  con- 
cept" [5],  it  is  possible  for  a  user  to  substitute  cartridge 
tape  for  disk  mass  storage  in  nearly  any  program,  this  prac- 
tice is  unsuitable  for  solution  of  large  systems  of  linear 
equations  out-of-core.   Previous  experience  with  less  demand- 
ing mass  storage  tasks  [4]  showed  degradation  of  system 
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performance,  on  large  problems,  where  tapes  were  required  to 
rewind  too  many  times.   Stretching  of  the  tape  initially 
introduced  parity  errors  and  later  caused  breakage  of  the 
tape  . 


Table  1:   Configuration  of  the  Hewlett-Packard 
Systems  45  Installation 


Required  Devices: 

Part  No.        Description 


HP  9845A 

HP  98032A 

HP  9  885M 

HP  9  885S 


Remarks 


Computer  with  CRT 

16  Bit  Interface  (9885A  Disk)  Option  85 

Flexible  Disk  Drive  (Master)  Drive  0 

Flexible  Disk  Drive  (Slave)  Drive  1 


Required  Options 
Numbe  r 

203 
310 
320 
330 
370 
500 


Descr ipt  ion 

64K  RAM 

Mass  Storage  ROM 

I/O  ROM  Left 

I/O  ROM  Right 

Graphics  ROM 

Internal  Thermal  Printer 


Remarks 


Select  Cede  13 
Select  Code  0 


Select  Code  Summary: 

Code  Device 


0 

3 

13 

16 


Printer 

Flexible  Disk  Drive  (Controller) 

Graphics  ROM 

CRT  Display 
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b.   Machine  Precision 

The  HP9845A  stores  floating  point  numbers  between 
9.99999999999  X  1099  and  1  X  10~"  in  absolute  value.   An  ex- 
tended calculation  range  of  9.99999999999  X  10511  to  1  X  10~511 
is  used  for  the  results  of  intermediate  calculations  [5]. 
Floating  point  numbers  are  described  as  full  precision  or 
short  precision.   Full  precision  floating  point  variables 
each  occupy  5  bytes  of  storage  and  are  used  exclusively  in 
all  HP  Enhanced  BASIC  Programs  presented  in  this  thesis. 

Chapter  3  of  Reference  [7]  discussed  the  estima- 
tion of  machine  roundoff  error  (as  part  of  a  discussion  on 
estimating  the  relative  error  in  the  solution  of  systems  of 
linear  equations  from  the  matrix  condition  number).   The 
worst  case  roundoff  error,  in  base  10  for  a  12  digit  machine 
is: 


io'12+1  =  0.00000000001 


Section  2.2  of  Reference  [7]  also  defines  a 
quantity  called  the  machine   epsilon  as  the  smallest  float- 
ing point  number  distinguishable  from  zero.   A  program, 
modeled  after  that  presented  in  the  reference,  found  the 
machine  epsilon  of  the  HP9845  to  be  7  .  2 75 9 5 7 6 14 15E-12 . 
c.   The  Enhanced  3ASIC  Language 

The  Enhanced  BASIC  Language  from  Hewlett-Packard 
has  all  the  features  of  the  most  powerful,  high  level  versions 
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of  BASIC  (Beginners  All-purpose  Symbolic  Instruction  Code) 
and  certain  extensions  designed  to  provide  optimum  use  of 
the  System  45.   References  [5]  and  [6]  contain  a  complete 
description  of  the  language  and  references  for  further  study 

As  with  most  versions  of  microcomputer  BASIC, 
the  language  is  completely  interpretive  in  nature.   That  is, 
no  compilation  or  conversion  of  source  code  to  machine  lan- 
guage is  accomplished  before  executing  the  program.   Each 
line  of  source  code  is  interpreted,  line-by-line,  as  execu- 
tion progresses.   Although  this  makes  it  possible  for  the 
program  to  have  such  incredible  flexibility  as  to  be  able, 
in  some  cases,  to  edit  itself  during  execution,  it  also  con- 
siderably slows  the  speed  of  execution.   On  the  positive 
side,  it  is  true  that  there  are  many  instructions  in  Hewlett- 
Packard  Enhanced  BASIC  which  can  reduce  the  amount  of  source 
code  required  over,  say,  the  same  program  in  FORTRAN  (e.g., 
the  MAT  series  of  instructions  for  matrix  manipulation  [5]). 

The  reader  who  is  familiar  with  FORTRAN  will 
have  little  trouble  following  the  programs  presented  in  this 
thesis.  The  table  below  lists  some  of  the  more  unusual  fea- 
tures of  the  language  and  should  be  enough  to  allow  complete 
understanding  of  program  flow.  For  further  information,  the 
reader  must  consult  the  references  given  above. 
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Table  2:   Selected  Members  of  the  Hewlett-Packard 
Enhanced  BASIC  Instruction  Set 


Def init  ions : 

" f ilenumber "  is  a  numeral  or  numeric  expression 
having  an  integer  value  between  1  and  10. 

"file  specifier"  is  a  character  string  enclosed  in 
quotation  marks  or  a  string  variable  which  is  a  valid 
and  complete  file  name  (including  mass  storage  unit 
specifier) . 

Stat  ement  s : 


ASSIGN//  f  ilenumber  TO  file  specifier 

The  ASSIGN  statement  places  the  filename  in  an 
internal  table  and  allows  the  file  to  be  referenced  by 
a  single  integer  number  in  PRINT//  or  READ//  statements 
[5]. 

BUFFER//  filenumber 

Sets  up  an  additional  256  byte  I/O  buffer  for  all 
READ'/  and  PRINT?/  statements  to  the  specified  filenumber, 
allowing  the  most  efficient  coupling  between  a  file  and 
its  associated  mass  storage  device  [6]. 

OPTION  BASE  1 

Specifies  that  the  lower  bound  of  all  array  variable 
subscripts  shall  be  1  (e.g.  Vector  (1)  is  allowed  but 
Vector  (0)  is  not)  [ 5  ]  . 

Y  y  y  =  F  N  X  x  x  x  x  x  (  a  a  ,  B  ,  C  c  c  ,  .  .  .  ) 

Call  to  a  function  subprogram  previously  defined  by 
a  DEF  FNXxxxxx  statement  [5].   In  FORTRAN  this  would  be 
YYY  =  XXXXXX  (AA,B,CCC, . . .)  where  XXXXXX  is  a  FUNCTION 
subprogram  or  an  intrinsic  function  set  up  with  a 
statement  function  definition. 

ZzzOOO : 

An  example  of  a  label  to  which  program  flow  may  be 
transferred  and  which  keeps  its  relative  position  within 
the  file  regardless  of  changes  in  line  numbering  [5]. 
Any  proper  statement  may  follow  the  colon. 


3EEP 


La  us 


es  a  tone  to  be  sounded  at  the  keyboard  [5]. 
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DISP  "Character  string" 

Causes  the  character  string  enclosed  in  quotation 
marks  to  be  printed  in  the  display  area  of  the  CRT  [5]. 

OVERLAP 

Allows  computation  and  I/O  to  occur  simultaneously. 
Most  beneficial  when  a  program  has  nearly  equal  amounts 
of  computation  and  I/O.   This  mode  is  disabled  by  the 
SERIAL  statement.   See  [5]  and  its  references  for  a 
discussion  of  the  OVERLAP  statements  impact  on  error 
traos  . 


When  attempting  to  reduce  the  run  time  of  a  pro- 
gram, it  is  important  to  recall  the  method  used  by  the  HP9845 
when  a  subroutine  or  function  call,  or  transfer  to  a  state- 
ment label  or  line  number,  is  encountered.   The  BASIC  inter- 
preter jumps  to  the  beginning  of  the  file  whenever  a  search 
for  one  of  these  items  is  instituted.   Run  speed  can  be  im- 
proved significantly  when  the  most  frequently  used  functions 
and  subroutines  are  placed  as  near  to  the  beginning  of  a  long 
program  as  possible.   Note  also  that  comment  lines  in  BASIC 
(statements  preceeded  by  REM  or  an  exclamation  point)  are  now 
a  liability  because  they  increase  the  length  of  a  program. 
The  amount  of  documentary  comments  should  be  limited  to  the 
minimum  necessary  to  understand  the  code. 

It  is  possible  to  overlay  subroutines,  by  using 
the  LINK  statement,  at  any  line  number  in  the  main  program 
(or  in  place  of  the  main  program  if  desired)  [5].   Overlaying 
reduces  the  amount  of  program  in  core,  allowing  the  storage 
of  more  variables  or  the  execution  of  programs  too  large  to 
fit  in  memory  in  a  single  piece. 
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d.   Disk  Mass  Storage  Considerations 

The  characteristics  of  the  disk  mass  storage 
hardware  used  with  any  microcomputer  have  a  significant 
impact  on  the  performance  of  the  system  in  finite  element 
calculations.   The  size  and  number  of  usable  physical  records 
on  the  disk,  size  of  the  file  directory  available  and  types 
of  file  structure  and  access  allowed  are  all  important  con- 
siderations in  attaining  optimum  use  of  out  of  core  storage. 

In  the  HP  9845  there  are  two  types  of  records  with 
which  the  programmer  must  be  concerned  [6].   A  phy  s  ical  record 
is  the  unit  of  storage  dealt  with  by  the  mass  storage  device 
itself  (and  is  usually  the  buffer  size  associated  with  the 
device).   Def  ined  records  are  the  smallest  addressable  units 
of  storage,  in  even  numbers  of  bytes,  which  can  be  individually 
accessed  by  the  user.   Whenever  possible  the  defined  record 
length  should  match  the  physical  record  length  of  the  device 
in  use.   In  the  HP  9845  this  matching  of  record  lengths  can 
provide  an  improvement  in  I/O  performance  of  two-  or  three- 
to-one  . 

Reference  [6]  describes  the  use  of  binary  DATA 
files  and  Direct  Memory  Access  (data  transfer  without  the 
use  of  buffers)  for  the  rapid  transfer  of  entire  arrays 
between  disk  mass  storage  and  core  memory.   This  method  can 
result  in  a  considerable  I/O  improvement  over  buffering  tech- 
niques, depending  on  the  storage  device  and  the  amount  of 
data  to  be  accomplished. 
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The  characteristics  of  the  HP  9885  Flexible  Disk 
Drive  are  listed  in  the  table  below. 


Table  3:   Selected  Characteristics  of  the 

Hewlett-Packard  9885  Flexible  Disk  Drive 

Storage  Capacity: 

Bytes  (total)-  ----------------  499,200 

Physical  records  (maximum)  ----------  1,950 

Physical  records  (useable)  ----------  1,901 

Files  (Director  size  limit)-  ---------  352 

Bytes  per  physical  record   ----------  256 

Accessing: 

Maximum  transfer  rate  (bytes  per  second)  -  -  -  46,000 


2 .   The  Apple-II  Plus  Personal  Computer 

The  Apple-II  Plus  Personal  Computer  is  manufactured 
by  Apple  Computer  Incorporated  of  Cupertino,  California  and 
is  designed  around  the  MOS  Technology  6502  microprocessor. 
The  maximum  addressing  range  is  2    (64K)  locations  on  a 
16  bit  parallel  bus.   The  data  bus   is  8  bits,  parallel  and 
bi-directional  and  the  CPU  clock  speed  is  1.203  MHz.   More 
technical  information  is  available  in  Reference  [8]. 
a.   System  Configuration 

The  baseline  Apple-II  Plus  central  processor 
does  not  vary  appreciably  from  installation  to  installation, 
although  enterprising  companies  have  recently  marketed  acces 
sories  which  can  actually  replace  the  6502  microprocessor 
with  the  more  popular  Z80  CPU  (and/or  increase  the  clock 
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speed)  .   The  possible  systems,  however,  and  the  languages 
available  on  those  systems  are  too  numerous  to  mention.   One 
special  purpose  and  six  general  purpose  "slots"  are  available 
at  the  rear  of  the  computer  for  connection  of  peripheral 
devices.   These  slots  are  supplied  power  and  have  access  to 
all  data  bus,  address  bus,  control  and  interrupt  lines. 

The  table  below  describes  the  configuration  of 
the  system  used  in  this  investigation. 


Table  4:   Configuration  of  the  Apple-II  Plus 
Personal  Computer  Installation 


Hardware : 

No.  used  Device  Remarks 

1  Apple-II  Plus  48K  ram 

1  Apple  Language  System  16K  RAM 

2  Apple  Disk-II  Interface  Card 

4  Apple  Disk-II  Floppy  Disk  Drives 

1  Apple  Parallel  Printer  Interface  Card 

1  Microline-80  Dot  Matrix  Graphics  Printer 

1  Sonv  Trinitron  Color  Video  Monitor 


So  f  twa  re : 

Apple  FORTRAN 
Apple  Pascal 


Slot  Connections: 

Slot  No.  Device  Remarks 

0  Language  System  Card 

1  Printer  Interface  Card 

Disk-II  Interface  Card       Volumes  11  and  12 
6    Disk-II  Interface  Card       Volumes  4  and  5 


b.   Machine  Precision 

The  Apple-II  Plus  stores  floating  point  numbers 
between  1.7  X  10     and  5.3  X  10"    in  absolute  value  and 
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integer  numbers  between  +32,768  and  -32,768  [8].   Floating 
point  variables  each  occupy  4  bytes  (two  machine  words)  of 
storage  and  integer  variables  occupy  2  bytes  (one  word)  . 
The  bit  allocation  for  a  floating  point  number  is: 
Bit:     31       30  ...  23       22  ...  16  and  15  ...  0 
Item:    sign       exponent  mantissa 

Floating  point  variables  having  a  maximum  of  6  decimal  digits 
may  be  stored. 

Using  the  methods  of  chapter  3  of  reference [ 7 ] ,  the 
estimate  of  worst  case  roundoff  error  is: 

10"6+1   =   0.00001 

The  machine  epsilon  (smallest  floating  poinc  num- 
ber distinguishable  from  zero)  was  investigated  using  the 
program  of  section  2.2  of  reference  [7]  and  was  found  to  be 
.119209E-06  . 

c.   The  Apple  FORTRAN  Language  and  the  Pseudo  Machine 

3 
References  [9]  and  [10]  discuss  the  Apple  FORTRAN 

language  and  explain  its  implementation  on  the  machine  and 
relation   to  the  Apple  Pascal  Operating  System.   Although 
it  is  not  necessary  to  know  any  Pascal  to  operate  the  system 
and  use  Apple  FORTRAN,  from  a  user  point  of  view  it  is  help- 
ful to  remember  that  this  version  of  FORTRAN  was  implemented 


3 
Apple  FORTRAN  is  a  registered  trademark  of  the  Apple 

Computer  Company. 
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exclusively  for  use  with  the  operating  system.   When  using  the 
Apple-II  described  in  this  thesis,  the  user  is  really  working 
on  a  virtual  or  "pseudo"  machine  which  has  been  created  by  the 
operating  system.   This  is  best  thought  of  as  a  software  gen- 
erated device  which  emulates  the  familiar  hardware  (e.g., 
registers)  found  in  a  real  computer. 

The  machine  language  of  this  virtual  microcomputer 
is  "pseudo  code"  (p-code)  which  resembles  a  true  machine  lan- 
guage.  P-code  is,  supposedly,  portable  among  all  computers 

which  operate  under  any  version  of  University  California,  San 

4 
Diego,  Pascal  (UCSD  Pascal).    This  would  also  require  that 

the  mass  storage  medium  upon  which  the  program  is  stored  and 

the  format  of  storage  be  compatible  with  both  machines  (e.g. 

an  Apple-II  5  inch  floppy  disk  is  not  a  compatible  storage 

medium  with  the  IBM  3033  at  the  Naval  Postgraduate  School). 

A  p-code  interpreter  program,  which  is  a  true  machine  language 

program  and  is  therefore  computer  dependent  converts  the  p- 

code  to  (in  the  Apple-II  case)  6502  instruction  codes.  All 

other  programs,  including  the  operating  system  itself,  are  in 

p-code  . 


UCSD  Pascal  is  a  registered  trademark  of  the  Regents  of  the 
University  of  California.   Use  thereof  in  conjunction  with  any 
goods  or  services  is  authorized  by  specific  licence  only  and 
is  an  indication  that  the  associated  product  or  service  has 
met  quality  assurance  standards  prescribed  by  the  University. 
Any  unauthorized  use  thereof  is  contrary  to  the  laws  of  the 
State  of  California. 


It  is  clear  then  how  the  FORTRAN  Language  option 
works  on  the  Apple-II  Plus.   The  Apple  FORTRAN  Compiler  is  a 
p-code  program  which  reads  TEXT  files  having  proper  FORTRAN 
syntax  and  in  turn  creates  from  them  a  pseudo  CODE  file  which 
is  compatible  with  the  virtual  machine.   It  is  also  clear 
that  Pascal  and  FORTRAN  routines  can  te  mixed,  through 
linking,  into  a  single  executable  CODE  file.   One  of  the 
original  advantages  of  this  microcomputer  system  was  that 
FORTRAN  was  available  and  that  the  great  volume  of  scien- 
tific program  which  are  currently  in  the  literature  were 
potentially  adaptable  to  the  Apple-II.   The  author  consid- 
ered that  the  need  to  use  Pascal  would  nullify  this  advan- 
tage.  Unless  the  reader  desires  to  use  the  Pascal  unit 
presented  in  Appendix  A  of  this  thesis,  the  Pascal  language 
need  not  be  mentioned  again. 

The  "proper"  FORTRAN  syntax  mentioned  above  is 
described  in  detail  in  Reference  [9].   Apple  FORTRAN  is 
based  on  the  American  National  Standard  Programming  Language 
FORTRAN  ANSI  X3.9-19785  (popularly  known  as  ANSI  subset 
FORTRAN  77) .   There  are  certain  deviations  from  the  stan- 
dard, both  deficiencies  and  extensions,  which  are  adequately 
dis  cus  sed  in  [ 9 ]  . 


Available  from  the  American  National  Standards  Institute, 
Inc.,  1430  Broadway,  New  York,  New  York,  10018. 
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Readers  familiar  with  FORTRAN  will  have  no  trouble 
following  the  programs  given  in  this  thesis.   Some  specialized 
statements,  called  compiler  directives,  will  be  unfamiliar  to 
most  users  but  are  discussed  in  detail  in  Appendix  A.   Input/ 
Output  statements  will  briefly  be  discussed  in  the  next  section. 

Double  precision  calculations  are  not  possible  in 
Apple  FORTRAN  as  it  is  presently  implemented, 
d.   Disk  Mass  Storage  Considerations 

Physical  records  on  the  Apple-II  diskette  are  called 
blocks  and  each  consist  of  512  bytes.   Page  25  of  Reference  [10] 
contains  a  discussion  of  the  technical  details  of  the  diskette 
and  its  storage  format  and  is  summarized  in  the  table  below: 


Table  5:   Selected  Characteristics  of  the  Apple 
Disk-II  Flexible  Disk  Drive 


Storage  Capacity: 
3ytes  (total)  ------------   143,360 

Tracks-  ---------------  35 

Sectors  ---------------  16 

Physical  records /bio cks  (maximum)  -  -  280 

Physical  records/blocks  (usable)-  -  -  274 

Files  (Directory  size  limit)-  -  -  -  -  77 


Disk  files  on  the  Apple-II  may  be  any  combination 
of  sequential  or  direct  access  and  formatted  or  unformatted. 
The  OPEN  statement  in  Apple  FORTRAN  corresponds  to  the  more 
familiar  DEFINE  FILE  statement.   The  format  of  the  OPEN 
stat  ement  is : 
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OPEN  (u,FILE=fname  [ , ACCESS= ' s t r ingl ' ] 
+  [ ,STATUS=,string2* ]  [ , FORM= ' s t r ing3 ' ] 

+  [,recl=rl]) 

, where , 

OPEN  (u 

is  required  and  must  appear  as  the  first 
argument  and  "u"  may  be  an  integer  variable  or 
expression. 

,  FILE" f name 

is  required  and  must  appear  as  the  second 
argument  and  fname  is  a  CHARACTER  expression 
of  the  complete  file  name  (including  disk,  name 
and  type  extension). 

.STATUS-' string! ' 

is  optional  and  'stringl'  must  be  one  of  the 
CHARACTER  constants  'OLD'  or  'NEW  (the  default 
is  'OLD* ) . 

.ACCESS-' string2 ' 

is  optional  and  'string2'  must  be  one  of  the 
CHARACTER  constants  'SEQUENTIAL'  or  'DIRECT' 
(the  default  is  'SEQUENTIAL'). 

,FORM»' string  3' 

is  optional  and  'string3'  must  be  either  the 
CHARACTER  constant  'FORMATTED'  or  'UNFORMATTED' 
(the  default  is  'FORMATTED'). 
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,RECL=rl 

required  and  only  allowed  for  direct  access 
files  where  "rl"  is  an  integer  variable  or 
exp  ress  ion . 

One  serious  limitation  is  present  when  using  direct 
access  files  in  Apple  FORTRAN  with  the  Pascal  Operating  System. 
When  locating  data  in  the  direct  access  diskette  file,  the 
system  apparently  uses  the  record  number  and  record  length 
given  by  the  user  to  compute  the  number  of  the  bytes  from  the 
beginning  of  the  file  to  the  data  to  be  read  or  written. 
Since  the  largest  integer  which  may  be  represented  in  the 
machine  is  32,768,  this  is  the  maximum  number  of  bytes  which 
may  be  present  in  a  direct  access  file.   For  single  direct 
access  files  the  maximum  number  of  records  is  (32763  /  RECL) , 
so  that  for  the  optimum  record  length  of  512  bytes,  a  direct 
access  file  mav  be  at  most  64  blocks  long. 
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II .   COMPARATIVE  SYSTEM  TESTING  USING  OUT  OF 
CORE  LINEAR  EQUATION  SOVLERS 


The  key  factor  in  determining  the  suitability  of  a  par- 
ticular microcomputer  system  for  solving  finite  element 
problems  is  the  machines  facility  for  solving  systems  of 
simultaneous,  linear,  algebraic  equations.   It  was  decided 
that,  before  proceeding  with  the  installation  of  a  code  on 
either  chosen  microcomputer,  a  standard  test  should  be  de- 
vised to  ensure  that  reasonalby  large  systems  of  equations 
could  be  solved.   As  a  point  of  interest,  it  was  also  de- 
sired to  compare  the  relative  speeds  of  the  two  machines, 
although  neither  would  be  excluded  from  the  study  solely 
because  of  low  solution  speeds.   Since  microcomputers  are 
so  much  less  expensive  than  larger  computers,  the  relative 
cost  of  CPU  time  is  also  proportionately  much  lower,  and  a 
slow  machine  which  is  sufficiently  accurate  and  has  the 
capability  of  handling  a  program  large  enough  for  a  variety 
of  element  types,  remains  acceptable.   3efore  discussing  the 
choice  of  solution  methods  used  in  the  testing,  it  is  help- 
ful to  recall  the  source  and  characteristics  of  the  systems 
of  equations  to  be  solved. 

A  familiar  and  most  common  method  of  finite  element 
analysis  is  the  displacement  based  method  (see,  for  example, 
[12])  which  gives  rise  to  a  system  of  constant  coefficient 
differential  equations  of  the  form: 
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MU  +  CU  +  KU  =  F 

, where , 

M  is  the  constant  coefficient,  mass  matrix  of  the 
finite  element  structure  and  usually  considers  the 
element  mass  as  a  consistent,  equivalent  masses  lumped 
at  appropriate  nodes. 

C  is  the  constant  coefficient,  damping  matrix  of  the 
finite  element  structure.   The  coefficients  of  this  matrix 
are  frequency  dependent  and  are  best  determined  from  the 
assembled  mass  and  stiffness  matrices  and  consideration 
of  experimental  results  [12]. 

K  is  the  constant  coefficient,  stiffness  matrix  of 
the  finite  element  structure.   This  matrix  is,  in  the 
direct  stiffness  method,  assembled  from  element  stiff- 
ness matrices,  which  in  turn  are  developed  from  strength 
of  materials  considerations  and  the  geometry  of  the 
chosen  element  type. 

U  is  a  matrix  of  time  dependent,  global,  nodal  dis- 
placements in  the  degrees  of  freedom  appropriate  for  the 
element  type  chosen  and  the  physical  structure  that  the 
element  assemblage  represents.  The  coefficients  in  U 
and  U  are  understood  to  be  the  global,  nodal  velocities 
and  accelerations,  respectively. 
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F  is  a  matrix  of  consistent,  applied  forces  acting 
in  the  global  degrees  of  freedom  and  arising  from  body, 
surface  and  concentrated  forces  applied  to  the  element 
structure  and  from  initial  stresses  (e.g.  thermal  or 
"lack-of-f it"  pre-stresses) .   In  the  most  general  case, 
these  forces  will  also  be  functions  of  time. 

For  static  analysis  M  and  C  effects  are  disregarded  and: 

K  U  =  F 

, where,  the  entries  of  U  and  F  are  independent  of  time  and 
the  system  of  differential  equations  reduces  to  a  system  of 
linear  algebraic  equations.   An  over  simplified,  but  helpful 
way  to  understand  this  system  is  to  recall   a  simple 
Newtonian  summation  of  forces,  including  D'Alembert's  dynamic 
forces,  on  a  single  mass  (i.e.,  F  =  ma);  the  two  situations 
are  analogous,  the  finite  element  system  being  more  general 
and  most  powerful.   Even  when  dynamic  effects  are  considered, 
direct  integration  methods  [12]  are  most  often  used  and 
require  only  the  solution  of  the  general  form  at  various 
instances  in  time,  such  that  one  is  again  solving  a  linear- 
ized algebraic  system.   The  equations  produced  by  the  above 
method  have  several  desirable  properties  [12-14]. 

Firstly,  the  system  of  equations  is  generally  well  con- 
ditioned and  positive  definite  [14].   For  this  class  of 
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problem  Gaussian  elimination  or  equivalent  techniques  (e.g. 

T 
see  LDL   factorization  in  [12])  are  the  most  efficient  and 

stable  [ 13]  . 

Secondly,  the  system  is  symmetric,  such  that  the  matrix 
of  coefficients  has  the  appearance  cf  having  been  reflected 
about  the  main  diagonal.   In  the  interests  of  reducing  the 
number  of  arithmetic  operations  and  the  storage  required  for 
solution,  only  the  main  diagonal  and  all  entries  above  (or 
below)  need  be  considered.   This  is  the  most  easily  pro- 
grammed improvement  ever  full  Gaussian  elimination  and 
requires  the  least  amount  of  additional  code  [13]. 

Thirdly,  the  coefficient  matrices  are  sparse  and 
tightly  banded,  if  proper  techniques  of  nodes  numbering 
are  used   [12,  14].    Sparsity  means  there  are  a  large 
number  of  zeros  in  the  off-diagonal  entries  of  the  matrix 
and  bandedness  means  that  most  non-zero  entries  are  con- 
centrated near  the  main  diagonal.   With  some  additional 
code,  the  solution  algorithm  can  be  designed  to  store  and 
work  only  with  elements  "below  the  skyline"  of  the  matrix 
[12],  and  to  have  some  reduction  in  storage  requirements  and 
solution  times.   Note,  however,  that  for  smaller  systems  of 
equations  and  for  microcomputers,  this  reduction  in  storage 
is  somewhat  offset  by  the  additional  code  required.   With 
an  interpretive  language  such  as  3ASIC,  some  of  the  speed 
advantage  is  also  offset  when  this  additional  code  must  be 
repeatedly  read. 
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Since  the  exactness  of  the  finite  element  solution  depends 
in  the  large  part  on  the  fineness  of  the  mesh  (greatest  pos- 
sible number  of  discrete  elements  to  model  the  continuous 
real  system),  it  is  not  uncommon  for  a  practical  problem  to 
have  from  hundreds  to  tens-of- thousands  of  elements  [12]. 
No  computer  presently  in  use  can  solve  the  largest  problems 
by  having  the  entire  system  of  equations  in  core  throughout 
the  entire  solution  phase.   Most  new,  general  purpose  codes 
work  with  out-of-core  solution  of  equations  and  store  the 
majority  of  the  system  on  magnetic  disk  or  tape  as  part  of  a 
more  complete  problem  database.  Any  useful  algorithm  must 
be  capable  of  solving  large  numbers  of  equations,  efficiently, 
out-of-core . 

A.   CHOICE  OF  EQUATION  STORAGE  SCHEME  AND  SOLUTION  ALGORITHM 

For  the  purposes  of  the  equation  solution  testing  phase 
of  this  thesis,  the  algorithm  of  Reference  [13]  was  selected. 
This  algorithm  was  perceived  to  be  well  suited  to  the  general 
test  role  for  three  reasons.   Firstly,  the  number  of  equations 
which  may  be  solved  is  theoretically  limited  only  by  the  amount 
of  out-of-core  storage  available  to  the  computer.   Secondly, 
the  amount  of  core  memory  required  for  solution  is  uncoupled 
from  the  number  of  equations  and  may  be  independently  selected 
when  choosing  the  "block"  size  (the  size  of  the  square  sub- 
matrices  into  which  the  complete  systems  is  divided).   No 
matter  how  large  a  percentage  of  the  microcomputer  memory  is 
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taken  up  by  the  solution  program,  a  sufficient  (if  not  optimum) 
amount  of  core  storage  should  remain  for  an  acceptable  block- 
size  to  be  chosen.   Lastly,  the  code  is  very  simple  to  follow 
and  requires  the  fewest  number  of  instructions  of  the  candi- 
date algorithms.   This  makes  a  program  easier  to  debug  and 
more  likely  to  fit  into  a  small  machine. 

3.   THE  BLOCK  SOLVER  METHOD  FOR  OUT  OF  CORE  SOLUTION  OF 
SYSTEMS  OF  LINEAR  EQUATIONS 

The  discussion  which  follows  parallels  that  of  Reference 
[13]. 

Consider  Figure  1,  which  shows  a  positive  definite,  sym- 
metric, sparse,  banded  matrix  in  an  associated  system  of 
linear  algebraic  equations.   The  terms  represented  by  "F" 
may  be  thought  of  as  loads  or  generalized  forces,  such  that 
"U"  and  "K"  might  then  represent  generalized  displacements 
and  equivalent  stiffnesses,  respectively.   The  general  e 1 e - 


F„"  and  "U  "  are  themselves  vectors  having,  say,  NS 


N 


components  and  the  general  terms  "K   "  are  NS-by-NS  square 
sub-matrices.   The  variable  NS  is  defined  as  the  block size, 
and  will  eventually  determine  the  amount  of  core  required 
to  solve  any  arbitrarily  large  system. 

Figure  2  illustrates  the  block  half -bandwidth  M  (or  MM 
or  Mm)  and  shows  that  portion  of  the  complete  stiffness 
matrix  which  must  be  stored,  for  this  algorithm,  to  effect 
a  solution.   Storage  for  fully  N-by-M  blocks  of  the  stiff- 
ness matrix  is  reserved,  even  though  this  is  more  than  is 
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absolutely  required  to  hold  one  complete  half-band.   The 
extra  blocks  are  used  for  storage  of  intermediate  results  of 
calculations  during  the  solution. 

The  program  requires  that  an  NS-by-NS  block  of  loads, 
each  containing  a  single  portion  of  the  complete  load  vector, 
be  stored  following  its  associated  row  of  stiffness  blocks. 
In  the  present  form  of  this  algorithm  only  one  "right-hand- 
side"  is  stored  in  each  load  block  and  most  of  the  space  in 
these  blocks  is  unused.   It  would  be  a  simple  matter,  once 
the  solution  procedure  is  understood,  for  the  reader  to  modify 
the  program  to  allow  up  to  NS  multiple  right  hand  sides  (and 
only  slightly  more  difficult  to  allow  a  greater  number) . 

It  is  recommended  that  all  blocks  (both  stiffness  and 
load)  be  stored  in  a  single,  large,  one  dimensional  array, 
unless  the  language  being  used  offers  special  advantages  for 
other  storage  schemes. 

The  block  solver  algorithm  itself  is  simply  a  form  of 
Gaussian  elimination  which  exploits  the  favorable  character- 
istics of  a  system  of  finite  element  equations.   The  reader 
who  is  interested  in  discovering  the  solution  scheme  for  him 
or  her  self  is  encouraged  to  complete  the  following  exercise, 
by  hand,  with  a  4-by-4  example: 

1.   Construct  the  first  equation  by  matrix  multiplication 
of  block  row  one  (K,  .,  -  K.-iM)  and  the  generalized  dis- 
placements (U,  -  U„)  . 
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Figure  2.   The  Banded  Character  of  the  Equivalent 
Stiffness  Matrix 
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2.   Solve  the  first  equation  for  U,  using  the  inver 

of  the  diagonal  block  of  the  current  row,  in  this 
.-1 


case  K 


11 


That  is 


Ul  +  Ku 


(F1  "  K12U2  -  K13D3 


-kimV 


Construct  the  next  equation  by  matrix  multiplication 
of  block  row  two  and  the  generalized  displacements. 
Substitute  the  expression  for  U..  ,  obtained  from  the 
first  equation,  into  the  second.   Collect  the  terms 
which  multiply  each  displacement  block  and  which 
involve  only  force  blocks. 

Use  the  collected  terms  to  reduce  load  and  stiffness 
blocks  as  shown  below: 


v       =    v 
"2     2 


T   -1 

"  K12K11 


K 
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=  K 


34 


T   -1 
-  K13K11K14 


or,  more  generally: 


* 

F   =  F 
M     M 


T   -1 
"  K1MK11F1 


K 


IJ 


-  ku  -  KiiKiiKu 
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The  expression  for  U1  obtained  in  step  2  must  be 

used  to  reduce  all  stiffness  coefficients  in  block 

rows  2  to  M  and  block  columns  2  to  M  and  all  load 

blocks  associated  with  those  block  rows. 

Once  the  reduction  in  step  6  is  complete,  the  second 

equation  is  solved  for  U       and  a  similar  reduction 

accomplished  on  block  rows  and  column  3  to  M+l  and 

their  associated  load  blocks. 

After  reduction  is  complete  through  the  Nth  block 

row,  the  Nth  equation  will  be: 


FN  "  KNN°N 


This  equation  is  solved  directly  for  U  ,  after  which 
back  substitution  is  used  to  solve  all  previous 
equations  in  the  usual  manner. 
Upon  completion  of  subroutine  SOLVE  the  original  stiff- 
ness matrix  has  been  destroyed  and  the  displacements  reside 
in  the  blocks  which  originally  contained  the  load  vector. 
The  algorithm  is  coded  using  five  subroutines  and  one 
function  subprogram  (or  statement  function  definition, 
depending  on  the  language  chosen)  .   Subroutine  SOLVE  is  the 
control  subprogram  which  carries  out  the  steps  outlined  in 
the  description  of  the  algorithm.   Subroutine  MULT  performs 
matrix  multiplication  of  two  arrays  which  have  been  stored 


in  unidimensional  form.   Subroutine  SYMINV  ensures  that  each 
diagonal  block  is  truely  symmetric,  attempts  to  detect  sin- 
gularity or  near  singularity  by  examining  the  diagonal  terms 
within  the  diagonal  block,  calculates  and  prints  a  condition 
number  and  inverts  the  block.   An  appropriate  warning  message 
is  printed  and,  if  necessary,  the  program  is  halted  whenever 
an  error  condition  exists.   The  condition  number  used  is  the 
product  of  the  maximum  column  sum  within  the  diagonal  block 
before  and  after  inversion  and  may  be  used  to  estimate  the 
solution  accuracy  [13].   Subroutines  RDDISK  and  WRDISK  are 
machine  dependent  subroutines  which  read  a  specific  record 
from  a  random  access  disk  file.   The  function  subprogram  used 
has  various  names  depending  on  the  version  (NTRACK,  FNTrk,...) 
but  is  simply  used  to  find  the  record  number  in  which  a  block 
from  block  row  I  and  block  column  J  is  stored.   Other  sub- 
routines required  to  generate  test  matrices  and  to  retrieve 
answers  and  control  programs  to  call  all  of  the  above  are 
discussed  in  the  next  section. 

1 .   Enhanced  3ASIC  Coding  on  the  HP-9845A 

The  block  solver  algorithm  was  coded  for  the  HP-9845A 
using  subroutines  SOLVE,  RDDISK,  WRDISK,  MULT,  SYMINV  and 
function  DEFinition  FNTrk.   Since  the  subroutines  are  not 
overlaid,  all  are  required  to  be  in  memory  at  the  same  time. 
All  transfers  within  the  controlling  subroutine  (SOLVE)  are 
made  to  labels  rather  than  statement  numbers,  so  that  the 
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program  may  be  RENumbered  at  the  convenience  of  the  user. 
In  order  to  decrease  the  run  time  of  the  program  REMarks 
were  kept  to  a  minimum,  a  BUFFER  was  used  and  the  most 
frequently  called  subroutines  were  placed  first  in  sequence. 
OVERLAP  mode  is  declared  upon  entry  into  subroutine  SOLVE 
to  allow  simultaneous  I/O  and  computation  and  SERIAL  mode 
is  declared  upon  exit  to  avoid  disabling  error  traps  in  the 
calling  program. 

The  HP-9845A  version  of  this  algorithm  is  given  in 
Appendix  B,  along  with  a  test  program  to  generate  various 
size  systems  of  suitable  equations  (subroutine  rest),  print 
the  answers   (subroutine  Answer)  and  to  mimic  the  FORTRAN 
EQUIVALENCE  statement  (subroutine  Equiv)  which  has  no  analog 
in  BASIC. 

A  sample  problem  and  program  output  is  provided  in 
Appendix  C. 

2 .   Apple  FORTRAN  Coding  on  the  Apple-II  Plus 

The  block  solver  algorithm  was  coded  for  the  Apple-II 
Plus  using  subroutines  SOLVE,  RDDISK,  WRDISK,  MULT,  SYMINV 
and,  within  SOLVE,  implicit  function  definition  NTRK(I,J). 
The  above  subroutines  are  all  compiled  into  one  ob j ec t / p-code 
file  by  imbeding  the  appropriate  SINCLUDE  statements  at  the 
end  of  SOLVE.   Dummy  dimensioning  (e.g.  AK1(1))  is  used  in 
all  subroutines  and  array  sizes  are  determined  in  the  calling 
program,  preferably  by  some  dynamic  means  using  a  single 
workspace . 
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In  order  to  handle  substantial  size  problems  in  the 
face  of  random  access  file  size  limitations,  subroutines 
RDDISK  and  WRDISK  are  set  up  to  handle  as  many  as  four  units 
(7-10),  each  of  which  must  have  been  pre-de fined /OPENed  in 
the  main  calling  program.   In  calling  RDDISK  or  WRDISK, 
NACTIV  is  the  desired  record  number,  C  is  the  array  to  be 
read  or  written,  NENTRY  is  the  one-dimensional  size  of  C 
and  NRECFY  is  the  number  of  records  per  random  access  file. 

The  Apple-II  Plus  version  of  this  algorithm  is  given 
in  Appendix  D,  along  with  a  series  of  subroutines  and  a  main 
program  to  OPEN  the  necessary  files  and  calculate  file  para- 
meters (subroutine  FILES),  generate  various  size  systems  of 
suitable  equations  (subroutine  TEST),  write  them  into  appro- 
priate files  (subroutine  DISKWR)  and  read  (subroutine  DISKRD) 
and  print  (subroutine  ANSWER)  the  answers.   The  control  pro- 
gram (program  THESIS)  is  also  included  in  Appendix  D  and  is 
organized  to  use  OVERLAYS  to  link  all  portions  of  the  test 
program  into  a  single  executable  CODE  file  and  to  manage  the 
workspace.   A  system  of  integer  pointers  (e.g.  NAK1 ,  NAK2 , 
...)  divide  up  a  workspace  of  3000  floating  point  variables 
according  to  the  chosen  block  size,  number  of  equations  and 
half-bandwidth . 

A  sample  problem  and  program  output  is  provided  in 
Append  ix  E . 
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3 .   Solut  ion  Times 

The  table  below  is  a  summary  of  all  the  test  runs 
made  on  each  machine.   Where  the  name  of  a  system  is  absent, 
that  particular  combination  of  parameters,  for  one  reason  or 
another,  was  not  attempted  (e.g.  the  blocksize  of  2  was  not 
tested  on  the  Apple-II  because  the  Disk  Drive  was  already 
continuously  busy  with  a  block  size  of  4) . 

For  comparison  purposes,  the  solution  times  obtained 
on  a  TEKTRONIX  4081  Graphics  System  Mini-computer  with  hard 
disk  drive  are  included  in  the  table.   For  more  information 
on  the  TEKTRONIX  version  of  the  block  solver,  the  reader  is 
referred  to  Reference  [15]. 

Sample  solutions  of  systems  of  160  equations  with 
true  half-bandwidths  of  64  appear  in  Appendices  C  and  E. 


Table  6:   Solution  Times  for  Standard  Equation 
Solving  T  rials 


Number  Block 

Block  Half- 

Rows  Bandwidth 

(NN)  (MM) 


16 
8 

4 

2 

1 

4 


16 
8 


S  qua  re  Recor  d 

Block  Length 
Size      In 

( N  S )  Bytes 


2 
4 

8 

16 

32 


32 

128 

64 

512 

256 

2048 

1024 

8192 

^096 

512 

256 


Computer 

Sy s  t  em 

Used 

HP9845A 

HP9845A 

App le-II 

HP9845A 

App le-II 

HP9845A 

Apple-II 

HP9845A 

Apple-II 

HP9845A 

App le-II 


Solut  ion 

Time 
(hr:min:sec) 

0 : 7 :45 . 44 
0:2:23.57 
0 : 7 :45  .13 
0: 3: 34  .93 
0:5:28.67 
0:1:53.10 
0:4:51.89 
0:2:41.67 
0:5:29.93 
0:0:52  .  33 
0:2 :50  .04 
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Table  6  (contd) 


Number 

Block 

Square 

Record 

Compute  r 

Solution 

Block 

Half- 

Block 

Length 

System 

Time 

Rows 

Bandwidth 

Size 

In 

Used 

(hr : min : sec) 

(NN) 
6 

(MM) 
3 

(NS) 
8 

Bytes 
512 

HP9845A 

0:2:22.92 

256 

Apple-I I 
TEK4081 

0:7:55.01 
0:0:09 

8 

4 

8 

512 

HP9845A 

0:5:05  .99 

256 

App le-II 
TEK4081 

0: 16:20.74 
0:0:23 

12 

5 

8 

512 

HP9845A 

0:11:48.64 

256 

Apple-II 
TEK4081 

0  :37:03.50 
0:0:40 

16 

6 

8 

512 

HP9845A 

0:22  :15  .91 

256 

Apple-II 
TEK4081 

1:08:37.68 
0:1:17 

20 

7 

8 

512 

HP9845A 

0:37  :14  .14 

256 

Apple-I I 
TEK4081 

1:51: 31.04 
0 :2 : 10 

20 

8 

8 

512 

HP9845A 

0:45 :47.61 

256 

Apple-II 
TEK4081 

2 : 13: 35 .89 

0:2:28 

C.   COMPARISON  AND  DISCUSSION  OF  RESULTS 

Problem  solution  times  for  both  the  HP-9845A  and  Apple-II 
Plus  greatly  exceeded  those  of  the  TEKTRONIX  4081  for  identical 
systems  of  equations.   Mote  that  the  TEKTRONIX  4081  never 
needed  2  1/2  minutes  for  any  problem  attempted  while  both 
of  the  test  systems  used  more  than  this  on  most  square  sets 
of  32  equations.   In  general  the  HP-9845A  needed  only  one 
third  the  solution  time  required  by  the  Apple-II  Plus. 

Although  it  was  originally  intended  that  the  TEKTRONIX 
4081  single  precision  solution  should  be  the  standard  for 
accuracy  (i.e.,  the  "exact"  answer),  it  was  discovered  that 
the  number  of  significant  digits  routinely  carried  in 
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HP-9845A  computations  exceeded  that  of  the  TEKTRONIX.   This 
more  precise  machine  word  resulted  in  the  answers   of  the 
HP-9845A  being  more  accurate  than  those  of  the  TEKTRONIX 
4081  doub le  precision  solution.   Using  the  HP-9845A  answers 
as  the  true  solution,  the  Apple-II  Plus  was  found  to  be  accu- 
rate to  three  significant  figures.   The  more  than  two  hours 
of  computation  (including  20  matrix  inversions)  required  of 
the  Apple-II  Plus  for  the  largest  problem  was  considered 
a  sufficiently  rigorous  test  of  machine  accuracy. 

The  first  five  runs  shown  in  Table  6  were  all  square 
systems  of  32  linear  algebraic  equations  in  32  unknowns. 
The  major  difference  in  these  problems  was  the  record  length, 
as  determined  by  the  square  blocksize  and  floating  point 
word  size,  although  the  problems  with  larger  blocksizes  were 
progressively  more  dense.   This  density  difference  was  negli- 
gible for  the  algorithm  chosen  because  the  skyline  within 
a  given  block  is  never  examined  nor  used  to  reduce  the  number 
of  arithmetic  operations.   It  was  expected  that  optimum  record 
lengths  would  be  those  which  matched  the  physical  record  size 
on  the  respective  disk  drive  system  and  so  result  in  the  fast- 
est Input/Output  performance.   In  seme  cases,  for  the  chosen 
algorithm,  this  particular  size  was  not  achievable  (i.e., 
when  the  disk  physical  record  length  in  bytes  divided  by  the 
floating  point  word  size  in  bytes  was  not  a  perfect  square). 
The  optimum  record  length  is  also  affected  by  the  relative 
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amounts  of  I/O  and  computation  in  a  program  and  the  ability 
of  a  machine  to  handle  both  simultaneously.   The  observed 
optimum  block  size  for  the  HP-9845A  and  Apple-II  Plus  was 
16  X  16.   The  remaining  test  combinations  were  run  to  gauge 
the  amount  of  time  to  be  expected  for  the  solution  of  larger 
systems  and  to  search  for  machine  peculiarities  which  might 
place  a  practical  limit  on  the  number  of  equations. 


43 


Ill .   IMPLEMENTATION  OF  REPRESENTATIVE  CODES 
ON  THE  TEST  SYSTEMS 

Upon  completion  of  initial  testing  of  the  equation  solv- 
ing capacity  of  the  HP-9845A  and  Apple-II  Plus  a  decision 
was  made  to  continue  the  study  by  implementing  a  finite  ele- 
ment code  on  these  systems.   The  plan  was  to  adapt  an  exist- 
ing code  to  solve  for  forces  and  displacements  in  space 
trusses  with  static  loading  and  to  which  other  element  and 
problem  types  could  later  be  added.   This  adapted  code  would 
be  re-programmed  and  modified  to  run  on  the  test  machines 
and  be  exercised  on  a  number  of  problems  of  suitable  diffi- 
culty.  It  was  considered  mandatory  that  the  source  program(s) 
use  an  out-of-core  equation  solver  and  desirable  that  some 
automatic  node  and  element  generation  features  be  provided, 
with  input  data  in  a  format  familiar  to  users  of  general 
codes.   Toward  this  end,  a  search  was  undertaken  to  find  an 
appropriate  source  program. 

A.        PREVIOUS    WORK    ON    THE    HP-9845A 

Reference  [3]  describes  Mallory's  success  in  modifying 
STAP  [12]  to  produce  SSAP-NPS  for  the  HP-9845A  with  tape 
mass  storage.   SSAP-NPS  solves  2-D  and  3-D  trusses  using  an 
in-core  equation  solver.   Although  some  thought  was  given 
to  the  possibility  of  modifying  it  to  solve  out-of-core  and 
to  adding  additional  elements,  time  constraints  finally 
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dictated  that  this  idea  be  abandoned  in  favor  of  developing 
a  code  for  the  Apple-II  Plus.   STAP-NPS  was,  however,  modi- 
fied by  changing  the  mass  storage  unit  specifiers  in  the 
program  to  select  the  newer  flexible  disk  drives  rather 
than  the  older  tape  units.   Test  runs  using  Mallory's  orig- 
inal problems  were  made  with  significant  reductions  in  solu- 
t  ion  time . 

3.   APPLE-II  PLUS  IMPLEMENTATION  OF  STAP-NPS 

The  challenge  of  implementing  a  finite  element  code  on 
the  Apple-II  Plus  really  involved  overcoming  the  limitation 
on  the  number  of  subroutine/function  calls  which  may  be 
made  by  a  program  (see  Appendix  A) .   The  random  access  file 
size  limit  (32768  bytes  maximum  length)  was  not  reached  in 
any  of  the  test  problems  run  in  for  the  thesis  (STAP-NPS 
will  presently  not  allow  the  stiffness  matrix  to  occupy 
multiple  files).   The  slow  speed  and  three  significant 
figure  accuracy  observed  in  the  Apple-II  were  already 
considered  in  the  first  phase  of  the  study. 
1 .   Software  Sources 

Many  of  the  popular  codes  in  use  today  are  exten- 
sions of  the  pioneering  work  of  Dr.  Edward  L.  Wilson  of  the 
University  of  California  at  3erkeley.   It  was  from  such 
sources  that  STAP-NPS  evolved. 
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The  main  program,  exclusive  of  element  dependent  sub- 
routines, was  taken  from  the  program  STAMOD   by  Nor,  while 
truss  element  subroutines  were  taken  from  Reference  [12]. 
2 .   Program  Development 

The  preliminary  development  of  STAP-NPS  for  the 
Apple-II  was  done  using  an  IBM  3033  mainframe  computer  and 
IBM  FORTRAN  at  the  W.  R.  Church  Computer  Center,  Naval  Post- 
graduate School,  Monterey,  California.   STAMOD  was  loaded 
from  an  existing  card  deck  while  subroutines  RUSS  and  TRUSS 
[12]  were  entered  by  hand.   Missmatch  problems  in  COMMON 
blocks  and  SUBROUTINE  calls  were  resolved  and  the  database 
structure  was  set  up.   The  resulting  program  was  then  con- 
verted from  double  to  single  precision  (the  Apple-II  has 
only  single  precision)  and  tested  on  several  problems. 

Up  to  this  point  STAP-NPS  was  one  single  program. 
The  next  step  was  to  split  the  program  up  into  separate  seg- 
ments, each  of  which  had  the  maximum  allowable  number  and 
depth  of  FUNCTION  and  SUBROUTINE  calls  and  which  would  fit 
in  the  Apple-II's  memory.   Estimates  of  what  would  "fit" 
were  based  on  previous  experience  and  careful  study  of  the 
number  of  variables  needed  in  memory  for  each  major  program 
phase.   It  was  finally  decided  to  break  the  program  into 
seven  parts . 


An  unpublished  program  used  in  supporting  studies  for 
"REALISATION  ET  ETUDE  D'ELEMENTS  DE  COQUE  DE  MINDLIN,"  Sopha 
Nor,  Doctoral  Dissertation,  Universite  de  Technologie  de 
Compiegne,  France,  29  June  1978. 
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3.   STAP-NPS 

The  FORTRAN  names  of  the  seven  PROGRAM  segments  are 
PROBLM,  ELEMS,  LOAD,  BLOCKS,  ASEMBL,  SOLVE  and  STRES .   Note, 
however,  that  these  programs  will  identify  themselves  by  a 
more  complete  title  when  interactively  communicating  with 
the  user  (e.g.,  ASSEMBLE  vice  ASEMBL  for  the  fifth  segment). 
Each  segment  has  certain  major  functions  which  it  is  expected 
to  accomplish  and  certain  subroutines  with  which  it  is 
associated. 

Program  PROBLM  determines  the  language  in  which  the 
output  is  to  be  printed  (both  French  and  English  are  avail- 
able), reads  control  parameters  and  nodal  point  information, 
automatically  generates  nodes  as  directed,  reads  and  stores 
boundary  conditions  for  global,  nodal  degrees  of  freedom  and 
determines  the  total  number  of  equations  to  be  solved. 

Program  ELEMS  calls  the  appropriate  element  subrou- 
tin  (only  a  truss  element  is  currently  available),  reads  and 
stores  connectivity  information,  automatically  generates 
elements  as  directed  and  calculates  the  characteristics  of 
the  global  stiffness  matrix. 

Program  LOAD  reads  in  loading  information  and  cal- 
culates and  stores  load  vectors. 

Program  BLOCK  considers  problem  size  and  user  recom- 
mended block,  size  and  determines  the  actual  blcoksize  to  be 
used.   Direct  access  record  length  and  actual  block  size  are 
then  used  to  develop  an  indexing  scheme  for  locating  any 
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coefficient  in  the  global  stiffness  matrix  within  its  random 
access  file.   The  number  of  columns  of  the  compact,  global 
stiffness  matrix  and  the  first  coupled  blocks,  per  block, 
are  calculated  and  stored. 

Program  ASEMBL  calls  the  appropriate  element  sub- 
routine and  calculates  and  stores  the  compact  form  of  the 
global  stiffness  matrix. 

Program  SOLVE  has  two  separate  modes  of  operation. 

T 
The  first  time  SOLVE  is  run  LDL   factorization  [12]  is  per- 
formed on  the  stiffness  matrix,  the  first  loadcase  vector  is 
reduced  and  multiplied  by  the  resulting  factored  form,  back 
substitution  is  carried  out  and  the  answers  (i.e.,  the  nodal 
displacements)  are  stored  for  calculation  of  stresses.   In 
subsequent  calls  successive  loadcase  vectors  are  retrieved 
and  used  to  calculate  their  respective  displacements. 

Program  STRES  calls  the  element  subroutine,  prints 
the  nodal  displacements  and  calculates  and  prints  the  ele- 
ven t  stresses. 

a.   Major  Variables 

The  major  variables  used  by  STAP-NPS  are  shown 
below : 
HUMNP        Number  of  nodal  points. 

Number  of  element  groups. 


NUMEG 


NUMI1AT  Number    of    materials. 


NUME 


Number  of  elements. 
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NDOF 
NEQ 
NLOAD 
NBLOCK 

ISTOTE 
ISTOH 


X(NUMNP) 
Y(NUMNP) 
Z(NUMNP) 
ID(NDOF*NUMN?) 


E(NUMMAT) 


AREA(NUMMAT) 


XYZ(6  *  NUME) 


Global  degrees  of  freedom  possible  per  node. 
Number  of  equations. 

Number  of  loads  within  a  given  loading  case. 
Number  of  blocks  into  which  the  global  stiff- 
ness matrix  is  divided  for  solution  purposes. 
Desired  number  of  coefficients  per  block. 
Actual  number  of  coefficients  per  block  as 
calculated  by  the  program  based  upon  the 
remaining  workspace  and  the  problem  size. 
X  coordinate  of  each  nodal  point. 
Y  coordinate  of  each  nodal  point. 
Z  coordinate  of  each  nodal  point. 
Contains  a  one  for  each  inactive  degree  of 
freedom  and  a  zero  for  each  active  degree  of 
freedom.   Zeros  are  later  changed  to  the 
equation  number  which  will  provide  the  solu- 
tion for  that  D.O.F.  and  node  and  the  ones 
are  changed  to  zeros  [12,  19]. 
Young's  modulus  for  each  material  group 
used  in  the  structure. 

X-sectional  area  for  each  material  group 
used  in  the  structure. 

Each  column  contains  the  cartesian  coordinates 
which  apply  to  a  single  element.   Rows  1-3 
and  4-6  contain  information  on  the  "i"  and 
"j"  ends  of  the  element,  respectively. 
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MHT(NEQ)         Vector  of  column  heights  in  the  global  stiff- 
ness matrix.   This  information  allows  a  com- 
pacted (having  fewer  zeros)  one-dimensional 
form  of  the  global  stiffness  matrix  to  be 
stored  and  accessed.   The  skyline  may  be  re- 
constructed from  the  MHT  vector  [12]. 

MAXA(NEQ  +  1)    Stores  the  addresses  of  the  diagonal  elements 

of  the  compacted  global  stiffness  matrix. 
MAXA(l)  is  always  one  and  MAXA(NEQ  +  1)  is 
always  the  total  number  of  elements  in  the 
compact  stiffness  matrix  plus  one. 

LM(6  *  NUME)     Each  column  of  this  array  holds  connectivity 

information  for  one  element.   Each  row  repre- 
sents a  local  element  degree  of  freedom  and 
each  entry  in  LM  is  a  global  equation  number 
taken  from  the  ID  array. 

MATP(NUME)        Contains  the  material  type  of  each  element. 

The  entries  in  this  vector  become  NUMMAT  for 
retrieving  the  Young's  modulus  and  x-sectional 
area  . 

R(NEQ)/V(NEQ)    A  single  Load  vector  describing  one  loadcase 

(i.e.,  a  "right-hand-side").   R  is  calculated 
from  the  individual  global,  nodal  loads  in 
FLOAD  and  the  ID  array. 
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NOD(NLOAD) 


IDIRN(NLOAD) 


ICOPL(NBLOCK) 


FLOAD(NLOAD)     Vector  of  individual,  concentrated  forces 

applied  to  global  nodes  within  one  loadcase. 
Vector  of  global  node  numbers  to  which  con- 
centrated forces  are  applied  within  one 
loadcase . 

Vector  of  principal  cartesian  coordinate 
directions  along  which  nodal  forces  are 
applied  within  one  loadcase. 
NCOLBV (NBLOCK)   Vector  containing  the  number  of  columns 

stored  within  a  given  block  of  the  com- 
pressed stiffness  matrix. 

Vector  containing  the  number  of  the  first 
block  which  is  coupled  to  another  particular 
block  of  the  compressed  stiffness  matrix. 
Various  blocks  of  the  compressed  stiffness 
matrix  during  the  assembly,  factorization 
and  back  substitution  phases. 

Vector  of  nodal  displacements  in  the  global 
degrees  of  freedom  (i.e.,  the  answer  to  the 
prob lem)  . 
b.   Management  of  Variable  Storage  Space 

In  programming  a  microcomputer  it  is  convenient 
to  use  one  single,  large  array  called  a  "workspace"  for 
storing  most  variables  which  will  be  READ  or  calculated  by 
the  program.   The  size  of  the  one  large  array  may  be  adjusted 
to  allow  the  maximum  number  of  variables  to  be  present  in 


AA  or  A  or  B 
(ISTOH) 


D(NEQ) 
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memory  along  with  the  program  instructions.   Pointers  are 
used  to  keep  track  of  and  locate  sub-arrays  within  the 
workspace  and  are  used  to  pass  the  locations  in  SUBROUTINE 
or  function  CALIs .   Because  the  Apple-II  does  not  use  the 
same  number  of  bytes  for  the  storage  and  representation  of 
integer  numbers  and  floating  point  numbers,  two  separate 
workspaces  (one  integer  and  one  real)  must  be  used.   Normal 
schemes  of  incrementing  the  pointers  differently  to  locate 
intermingled  floating  point  and  integer  sub-arrays  cannot 
be  used . 

The  integer  workspace  is  the  array  IA  which  is 
dimensioned  MTOTI  and  is  contained  in  COMMON  block  IWORK. 
The  real  workspace  is  the  array  A  which  is  dimensioned  MTOTR 
and  is  contained  in  COMMON  block  RWORK .   Integer  and  real 
pointers  begin  with  NI  and  NR,  respectively,  and  change 
their  true  position  within  the  workspaces,  during  each  major 
phase  of  the  problem,  as  the  variables  they  represent  are 
shuffled  about.   When  the  workspaces  are  passed  to  subrou- 
tine RWCOMN  the  pointers  NENDI  and  NENDR  mark  the  highest 
numbered  elements  of  IA  and  A  which  must  be  written  into 
and,  in  the  next  segment,  read  from  mass  storage.   The  fol- 
lowing two  tables  show  the  pointer  allocations  throughout 
the  solution  process: 
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The  pointers  are  modified  in  several  programs  and 
subroutines  in  STAP-NPS .   Each  time  the  numerical  value  of  the 
highest  active  pointer  is  increased  a  check  must  be  made  to 
ensure  that  no  element  of  the  workspace  will  be  stored  outside 
its  dimensioned  area.   The  magnitude  of  the  topmost  active 
pointer  is  compared  to  MTOTI  or  MTOTR,  as  appropriate,  and 
if  the  required  value  is  too  large  subroutine  ERROR  is  called. 
Subroutine  ERROR  prints  a  diagnostic  message  showing  the 
amount  by  which  the  offending  pointer  exceeds  the  memory 
available  and  the  problem  phase  which  cannot  be  accomplished. 
The  user  then  has  the  option  of  attempting  to  set  MTOTI  and 
MTOTR  closer  to  the  machine  limits,  rearranging  the  mesh  into 
smaller  pieces  (e.g.,  by  using  multiple  element  groups)  or 
abandoning  the  Apple-II  for  a  larger  machine  and  program. 

The  first  six  segments  of  STAP-NPS  all  have 
entry  and  exit  diagnostics  which  will  display  the  contents 
of  the  integer  and  floating  point  workspaces.   The  printing 
of  these  diagnostics  may  be  declined  by  the  user.   They  were 
first  included  in  the  programs  to  assist  in  debugging  and 
were  allowed  to  remain  so  that  users  could  more  easily 
modify  or  add  to  the  code.   These  diagnostics  could  also 
be  of  value  to  the  beginning  student  of  the  finite  element 
met  hod  . 
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c.   The  Database  and  Communication  Between  Program 
Segmen  ts 

Communication  between  programs  was  accomplished 
by  including  a  new  subroutine  (RWCOMN)  in  each  segment. 
This  subroutine  writes  the  contents  of  the  entire  COMMON 
block  (which  is  the  same  group  of  variables  in  each  segment) 
and  all  other  necessary  integer  and  real  variables  into  an 
unformatted,  serial  disk  file.   All  segments,  except  the 
first,  resume  solving  the  problem  after  calling  RWCOMN  with 
an  argument  of  "1",  thereby  receiving  the  necessary  data 
from  previous  calculations  and  READs .   All  segments  call 
RWCOMN  with  an  argument  of  "2"  as  their  final  action  before 
closing  file  ICOM  and  stopping.   The  old  ICOM  file  from  the 
previous  segment  is  deleted  just  before  the  new  one  is 
written,  well  after  any  run  time  errors  could  affect  its 
contents . 

Program  STAMOD  was  already  organized  around  six 
FORTRAN  I/O  units  which  naturally  evolved  into  the  remaining 
database  files  for  STAP-NPS.   These  database  files  are: 
UN 
This  sequential,  formatted  file  contains  the  input  data 
for  the  problem  to  be  solved.   UN  is  an  internal  name  only; 
the  actual  name  of  this  file  is  entered  interactively  by 
the  user.   (If  CONSOLE:$  is  entered  then  UN  becomes  the 
keyboard.)   It  is  created  by  the  user  with  the  Apple  E(ditor 
before  program  PROBLM  is  run  and  is  read  by  programs  PROBLM, 
ZLEMS  and  LOAD. 
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IOUT 

This  sequential,  formatted  file  contains  all  printed 
output  for  the  problem  to  be  solved.   IOUT  is  an  internal 
name  only;  the  actual  name  of  this  file  is  entered  inter- 
actively by  the  user.   (If  a  file  other  than  PRINTER:$  or 
CONSOLE:$  is  specified,  the  contents  must  be  examined 
before  the  next  segment  is  run  or  be  lost.) 
IDTAP 

This  sequential,  unformatted  file  contains  the  ID  array. 
It  is  created  in  program  PROBLM  and  read  in  program  SOLVE. 
IELMNT 

This  sequential,  unformatted  file  contains  an  indexing 
parameter  (MIDEST) ,  the  element  type  (NPAR1)  and  NUME , 
NUMMAT,  E,  AREA,  XYZ ,  LM ,  and  MATP  for  each  element  group. 
It  is  created  in  program  ELEMS  and  read  in  programs  ASEMBL, 
SOLVE  and  STRESS. 
ILOAD 

This  sequential,  unformatted  file  contains  R  vectors 
from  the  various  load  cases.   It  is  created  in  program 
LOAD  and  read  in  program  SOLVE. 
IRIG 

This  random  access,  unformatted  file  contains  the  con- 
densed, and  later  factored,  form  of  the  global  stiffness 
matrix.   It  is  created  in  program  ASEMBL  and  read  in  program 
SOLVE.   The  record  length  is  user  selectable  by  changing 
the  variable  LONG  in  program  PROBLEM. 
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ICOM 
This  sequential,  unformatted  file  contains  the  active 
portions  of  the  real  and  integer  workspaces  and  the  contents 
of  all  other  COMMON  blocks.   It  is  created  or  re-written 
in  each  program  segment  and  all  but  the  first  read  this 
file  through  subroutine  RWCOMN. 
4 .   Running  STAP-NPS 

In  order  to  "X(ecute"  the  programs  in  STAP-NPS  the 
following  disks,  at  a  minimum,  are  needed: 

F0RT1 : 

Contains  APPLE  FORTRAN  (Pascal)  Operating  system. 
Also  contains  the  F(iler  and  E(ditor  for  creating  input 
data  files  for  the  problem  to  be  solved. 

DATA1:  (or  DATA3 : ) 

Contains  executable  CODE  files  for  the  first  five 
segments /p rograms  in  STAP-NPS.   Also  contains  the 
input  data  TEXT  files  for  the  four  problems  used  for 
original  program  testing. 

DATA2:  (or  DATA4:) 

Contains  executable  CODE  files  for  the  last  two 
segmen ts /p ro grams  in  STAP-NPS. 

TWO  SCRATCH  DISKS 

Any  two  formatted  disks  having  sufficient  space 
for  database  files.   These  two  disks  must  NOT  be  write 
protected.   Normal  program  execution  will  not  damage 
pre-existing  files. 
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To  actually  run  the  programs: 

(1)  Load  volume  4  with  F0RT1:,  volume  5  with  DATAl: 
and  volumes  11  and  12  with  scratch  diskettes. 

(2)  Turn  on  the  monitor  and  other  peripherals. 
Energize  the  Apple-II  Plus  and  observe  initial 
program  loading  of  the  FORTRAN  Operating  System. 

(3)  Enter  the  F(iler  and  set  the  date.   Enter  the 
E(ditor  and  create  file  UN  using  the  format 
shown  in  the  following  subsection. 

(4)  X(ecute  programs  #5 : ONE . CODE  through  #5 : F IVE . CODE . 
These  are  simply  the  executable  p-code  files  of 
programs  PROBLM  through  ASEMBL.   Interact  and 
provide  information  as  requested. 

(5)  Remove  DATAl:  from  volume  5  and  replace  it  with 
DATA2 : 

(6)  X(ecute  programs  V 5 : S IX . CODE  and  # 5 : SEVN . CODE 
(executable  p-code  files  of  SOLVE  and  STRESS) 
alternately  until  all  loadcases  have  been  solved 
and  the  displacements  and  stresses  printed. 

a.   Input  Data  Card  Formats 

Input  data  cards  are  actually  lines  of  data  in 
database  file  UN.  The  user  is  required  to  construct  this 
file  using  the  E(ditor,  prior  to  running  STAP-NPS. 
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Language  Option  Card  (20A4) 

Note      Columns      Variable       Entry 
(1)        1-4  Talk  FRAN  for  French 

ENGL  for  English 
Notes: 

(1)  This  card  must  be  included.   Default  language  is 
Eng lish . 
Heading  Card  (20A4) 

Note      Columns      Variable       Entry 
(1)         1-30  HED  Problem  title 

Notes : 

(1)  This  card  must  be  included. 
Master  Control  Card  (115,611,114,315) 

Note      Columns      Variable       Entry 

1-5  NUMNP  Total  number  of 

nodal  po  int s 
(1)         6-11  IDOF  EQ.  0  free 

EQ.  1  fixed 
12-15  NUMEG         Number  of  element 

groups 
16-20  NLCASE        Number  of  load  cases 

21-25  MODEX         EQ .  0  Data  check 

EQ  .  1  Solut  ion 
26-30  ISTOTE        Desired  number  of 

coefficients  per 
block 


Notes : 

(1)  For  two  dimensional  trusses  enter  "001111." 

For  three  dimensional  trusses  enter  "000111." 
Default  is  "000000"  which  will  use  excessive 
diskette  storage  and  memory  for  the  solution. 
Nodal  Point  Definition  Cards  ( Al , 14 , Al , 14 , 515 , 3F10 . 0 , 15) 
Note      Columns      Variable       Entry 

1  IT  Blank  for  X,Y,Z, 

EQ .  X  Cylindrical 

(1)  2-5  N  Node  number 

(2)  6  JPR  Print  control 

see  Note  (2) 

(3)  7-10  IDT(l)         EQ.  0  free 

EQ.  1  fixed 
(EQ.-l  fixed  for 
automatic  gener.) 

(3)  11-35  IDT(2)-       EQ.  0  free 

IDT(6)         EQ  .  1  fixed 

(EQ.-l  fixed  for 
automatic  gener.) 
35-45  X  X  coordinate  of 

node  N 

(4)  46-55  Y  Y  coordinate  or 

radius  -  node  N 
(4)       56-65  z  Z  coordinate  or 

theta  -  node  N 
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Note      Columns      Variable       Entry 
(5)        66-70        KN  Node  increment 

for  automat  ic 
generat  ion 
Notes : 

(1)  Need  not  be  input  in  numerical  order  but  eventually 
must  read  in  nodes  1  through  NUMNP .  Last  node  read 
in  must  be  node  NUMNP. 

(2)  Print  control  character  only  read  from  card  for 
node  number  1.   Codes  are: 

Blank  -  no  print  supression 
EQ.  A  -  supress  list  of  node  coordinates 
EQ .  3  -  supress  list  of  equation  numbers 
EQ.  C  -  supress  both  of  the  above 

(3)  These  are  local  degrees  of  freedom  and  should  be 
fixed  where  global  conditions  so  require.   When  two 
consecutive  cards  are  being  used  to  automatically 
generate  elements,  all  fixed  degrees  of  freedom 
must  be  set  to  "-1"  vice  "1". 

(4)  Cylindrical  coordinates  only  if  IT  is  non-blank. 

(5)  Nodes  generated  in  the  sequence: 
N0DE1,N0DE1+(1*KN1)  ,  NODE  1+ ( 2 *KN1)  ,  .  .  .,  N0DE2 
where  N0DE1  is  N  on  the  first,  KN1  is  KN  on  the 
first  and  N0DE2  is  N  on  the  second  of  two  consecu- 
tive Nodal  Point  Definition  Cards.   Default  value 
of  KN  is  "I". 
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Element  Data  Block  Marker  Card  (1A4) 

Note       Columns      Variable       Entry 
(1)         1-4  BLOCK         ELEM 

Notes  : 

(1)  This  card  must  be  included.   It  is  the  marker 

for  Program  ELEMS  to  begin  reading  element  data 
from  file  UN  after  unit  rewind. 
Material  Group  Specification  Cards  (I5,2F10.0) 
Note      Columns      Variable       Entry 

1-5  N  Material  group 

6-15  E  Young's  Modulus 

16-25  AREA  Truss  element 

x-sectional  area 
Notes : 

None  . 
Element  Data  Cards  (515) 

Note      Columns      Variable       Entry 
(1)         1-5  M  Element  number 

6-10  II  Node  number  at 

one  end 

Node  number  at 
other  end 
Material  type 
Increment  for 
automatic  element 
generat  ion 


(2) 


11-15 

16-20 

21-25 


J  J 

MTYP 
KG 
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Notes : 

(1)  Put  in  ascending  order  beginning  with  number  "1". 

(2)  Missing  elements  are  filled  in  by  using  the  material 
type  of  the  last  card  before  the  missed  element  and 
incrementing  the  element  numbers  by  one  and  the 
node  numbers  by  KG  until  the  gap  is  closed. 

Loading  Data  Block  Marker  Card  (1A4) 

Note      Columns      Variable       Entry 
(1)         1-4  BLOCK  LOAD 

Notes : 

(1)  This  card  must  be  included.   It  is  the  marker  for 
Program  LOAD  to  begin  reading  loading  case  data 
from  file  UN  after  unit  rewind. 
Load  Case  Control  Cards  (215) 

Note      Columns      Variable       Entry 

(I)         1-5  LL  Loading  case  number 

6-10  NLOAD         Number  of  loads 

within  the  load 
Notes : 

(1)  Enter  is  ascending  order  beginning  with  "1". 
Concentrated  Load  Cards  (2I5,F10.0) 

Note      Columns      Variable       Entry 

1-5  NOD  Node  to  which  load 

is  applied 
(1)         6-10  IDIRN  Pinciple  direction 

11-20  FLOAD         Load  magnitude 
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Notes : 

(I)  Positive  directions  only: 

EQ.  1   -   +X 

EQ  .  2   -   +Y 

EQ.  3   -   +Z 
For  negative  directions  change  the  sign  of  FLOAD. 
Language  Option  Cards,  Heading  Cards,  Master  Control  Cards 
and  Data  Block  Marker  Cards  will  appear  only  once  in  every 
file.   Other  cards  will  appear  some  number  of  times,  depen- 
dent upon  the  particular  problem  being  solved. 


C.   TESTING  OF  STAP-NPS 

STAP-NPS  was  tested  on  four  problems,  the  formatted 
input  files  of  which  appear  in  Appendix  G.   Problems  were 
chosen  for  their  varying  degree  of  difficulty,  the  combina- 
tion of  features  which  they  tested  and  the  availability 
of  a  previous,  reliable  solution. 

Truss  number  one  was  a  simple  three  member  problem, 
easily  solvable  by  hand,  which  was  used  to  debug  the  program. 
Truss  number  two  had  more  nodes  and  elements  and  tested  the 
multiple  load  case  feature.   Truss  number  three  was  a  large 
crane  with  multiple  materials  and  automatic  node  generation. 
Truss  number  four  was  a  problem  taken  from  Mallory's  thesis 
[3]  and  tested  the  automatic  element  generation  feature  of 
the  program  and  three  dimensional  operation.   The  cylindrical 
coordinate  option  and  print  supression  features  have  never 
been  tested. 
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The  English  language  program  output  from  solving  truss 
number  four  is  given  in  Appendix  H.   The  output  of  the  first 
segment  for  the  same  problem  with  French  as  the  chosen  lan- 
guage has  also  been  included. 
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IV.   CONCLUSIONS 

In  concluding,  two  separate  questions  will  be  addressed. 
Firstly,  what  general  techniques  should  be  used  when  tailor- 
ing a  finite  element  code  to  a  microcomputer?  Secondly,  how 
well  suited  to  finite  element  work  were  the  two  test  systems 
and  what  limitations  were  encountered  in  each? 

When  using  a  microcomputer  system,  programs  should  be 
designed  around  a  complete  data  base.   Files  should  contain, 
in  so  far  as  possible,  only  one  type  of  information  and  have 
a  format  general  enough  for  a  variety  of  purposes  (e.g., 
nodal  point  information  for  use  by  both  stress  calculation 
and  graphics  mesh  plotting  subroutines).   Single  workspaces 
with  dynamic  dimensioning  (i.e.,  pointers  that  change  with 
the  problem)  are  recommended  for  the  ease  with  which  avail- 
able memory  may  be  managed.   Overlay  schemes  for  subroutines 
should  be  used  to  ensure  that  a  minimum  amount  of  memory  is 
taken  up  by  program  instructions.   In  choosing  the  extent  of 
overlaying  to  be  used,  study  is  required  to  determine  the 
optimum  balance  between  I/O  and  computation  and  the  amount 
of  core  needed  for  variable  storage. 

Never  select  microcomputer  operating  system  software 
(regardless  of  the  language  in  use)  which  has  low  limits  on 
the  number  of  user  subroutines  accessible  to  main  programs. 
Subroutines  or  functions  should  be  able  to  be  individually 
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overlayed  even  though  they  may  reside  in  large  libraries  of 
many  such  routines.   Beware  of  eight  bit  microcomputers 
whose  small  word  size  leads  to  accuracy  problems  (when 
double  precision  software  is  not  available)  and  limits  mass 
storage  file  size  (i.e.,  where  the  operating  system  chooses 
to  keep  track  of  total  bytes  per  file). 

The  Hewlett-Packard  System  45  Desktop  Computer  exhibited 
only  three  negative  attributes  during  this  study.   The  first 
of  these  was  its  slow  speed  (a  problem  which  has  no  practical 
solution).   The  second  was  the  lack  of  a  FORTRAN  compiler. 
The  third  was  the  relatively  high  cost  of  the  HP-9845A  com- 
pared to  other  micro-  and  even  mini-  computers.   Although 
the  System  45  is  extremely  accurate  and  capable,  it  is  much 
in  demand  for  graphics  related  work  (for  which  it  is  perhaps 
best  suited)  and  it  is  probably  not  cost  effective  to  devote 
so  expensive  a  machine  to  hours  of  slow  and  deliberate  num- 
ber crunching.   Organizations  who  wished  to  use  the  HP-9845A 
for  finite  element  calculations  would,  however,  experience 
few  difficulties. 

The  Apple-II  Plus  Personal  Computer,  with  the  operating 
system  and  hardware  configuration  previously  discussed,  is 
not  a  suitable  tool  for  serious  finite  element  work.   The 
reader  will  recall  that: 
(1)  The  system  is  too  slow,  taking  over  two  hours  to  solve 

160  equations  having  a  half -bandwidth  of  64. 
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(2)  There  is  a  limit  of  seven  compiler  calls  ($USES  statements) 
to  external  units  for  other  subroutines  or  functions. 

(3)  Random  access  files  are  limited  to  32,768  bytes  in 
total  length . 

(4)  The  machine  word  is  too  small  and  double  precision  is  not 
available.   This  causes  a  loss  of  accuracy  in  practical 
calculat  ions . 

(5)  Using  more  than  one  implicit  function  definition  per 
program  causes  multiple  subroutines  with  the  name  DUMMY 
to  be  generated  and  results  in  compile  time  errors. 

(6)  Many  machine  features  cannot  be  used  with  the  FORTRAN 
language  alone. 

If  the  Apple-II  must  be  employed  it  is  recommended  that 
another  operating  system  and  central  processor  chip  be  used. 

The  author  looks  forward  in  anticipation  to  the  wedding 
of  the  finite  element  method  and  the  microcomputers  of  the 
future  and  to  the  powerful  vehicle  of  discovery  which  will 
then  be  widelv  available. 
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APPENDIX  A 
SEPARATE  COMPILATION  OF  APPLE  FORTRAN  PROGRAM  SEGMENTS 

Due  to  the  small  size  of  the  Apple-II  Plus  memory,  it 
is  often  necessary  to  compile  larger  programs  in  separate 
pieces  and/or  to  overlay  certain  subroutines  during  execu- 
tion of  a  large  main  program.   This  is  made  possible  in 
Apple  FORTRAN  through  the  use  of  compiler  directives  [9]. 

Compiler  directives  are  used  to  communicate  information 
to  the  compiler  via  the  FORTRAN  TEXT  file.   These  statements 
all  begin  with  a  dollar  sign  in  column  one  and  some  have 
special  requirements  as  to  their  location  within  a  file 
(e.g.  3USES  statements  must  be  placed  at  the  beginning  of  a 
TEXT  file  and  may  only  be  preceeded  by  Comment  1 ines )  .   For 
details  about  the  majority  of  compiler  directives  statements 
see  Reference  [9].   A  brief  commentary  on  some  of  the  state- 
ments appears  in  the  table  below: 

Table  9:   Compiler  Directives  in  Apple  FORTRAN 


SINCLUDE  filename 

Tells  the  compiler  when  this  statement  is 
encountered  to  immediately  compile  the 
contents  of  the  FORTRAN  TEXT  file  filename. 
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Table  9  (contd) 

$XREF 

This  statement  produces  a  cross-reference 
listing  of  the  compilation.   Note  that  the 
extra  information  generated  by  this  statement 
reduces  the  size  of  programs  which  may  be 
compiled.   If  not  absolutely  required  this 
directive  should  be  avoided. 

$EXT  [type]  FUNCTION  name  Iparais 
$EXT  SUBROUTINE  name  rfparams 

Used  to  identify  assembly  language   program. 

calls  and  not  used  in  this  thesis. 


The  format  of  the  $USES  control  statement  is : 

$USES  unitname  [  IN  filename  ]  [  OVERLAY  ] 

where  brackets  enclose  optional  items  and  capitalized  items, 
when  appearing,  must  be  exactly  as  shown. 

The  character  string  "unitname"  is  a  Pascal  term  which 
actually  represents  a  collection  of  P-code  "procedures"  or 
"functions"  (FORTRAN  subroutines  or  FORTRAN  function  sub- 
programs previously  compiled  on  the  Apple)  with  their 
associated  compiler  directives.   This  collection  is  com- 
piled as  a  single  entity  and  remains  together  during  any 
OVERLAY.   The  great  disadvantage  in  this  technique  is 
that  more  code  will  be  called  into  core  than  is  actually 
needed  whenever  unused  subroutines  reside  in  the  same  unit 
as  one  which  is  desired.   If  one  decides  to  simply  put  a 
single  subroutine  in  each  unit  to  circumvent  this  problem, 
there  remains  a  limit  on  the  maximum  number  of  $USES  state- 
ments which  are  allowed  by  the  compiler  in  any  one  program. 
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The  character  string  "filename"  includes  all  of  the 
information  necessary  for  the  compiler  and  linker  to  find 
the  file  and  usually  includes  the  diskette  name  or  volume 
number  followed  by  a  colon,  the  file  name  followed  by  a  per 
iod  and  the  extension  CODE.   If  the  filename  is  not  included 
the  compiler  and  linker  will  expect  to  find  the  unit  in  the 
file  7/4:  SYSTEM.  LIBRARY  . 

The  character  string  OVERLAY  tells  the  compiler  and 
1  inker  that  the  unit  named  is  to  be  loaded  into  core  only 
when  called  and  is  to  remain  in  memory  only  while  the  pro- 
cedure of  interest  is  in  use. 

A  maximum  of  seven  $USES  statements  may  appear  in  any 
one  program  or  subprogram  being  compile.   Whenever  more  than 
this  number  appear  in  a  file,  compile  time  error  number  205 
is  generated  and  the  system  usually  "hangs"  and  must  be 
re-booted.   This  is  a  real  inconvenience  when  setting  up  a 
program  in  modular  form  and  assigning  all  repetitive  tasks 
to  subroutine  or  in  setting  up  a  large  overlay  scheme.   For 
example,  of  one  uses  only  one  level  of  depth  in  program 
development  (e.g.  the  main  program  calls  all  subroutines 
and  no  subroutine  calls  any  other)  then  the  program  may  be 
broken  into  only  seven  parts,  some  of  which  may  not  be  small 
enough  to  allow  storage  for  the  desired  number  of  variables 
to  reside  in  core.   If  graphics  are  required  in  any  sub- 
routine, one  of  the  seven  SUSES  statements  will  probably 
have  been  used  to  link  the  TURTLEGRAPHICS  unit  from  the 


#4: SYSTEM. LIBRARY .   For  programs  where  subroutines  call  or 
overlay  other  subroutines  the  situation  is  much  more  complex 
and  is  best  illustrated  by  the  example  in  the  following 
section . 

Reference  [9]  discusses  the  use  of  the  $USES  compiler 
directive  when  more  than  one  level  of  depth  is  used  in  pro- 
gram development,  but  gives  no  actual  examples.   Since  the 
solution  of  this  problem  is  so  critical  to  the  success  of 
implementing  a  large  program  on  the  Apple-II  Plus,  it  was 
decided  to  examine  a  complex  case  of  subroutine  inter- 
connec  t  ion . 

The  tables  which  follow  contain  actual  files  which  were 
successfully  compiled,  linked  and  tested  on  the  Apple-II. 
Although  the  files  are  real,  they  are  simply  presented  as  a 
hypothetical  case  which  exercises  the  compiler  and  linker  to 
the  maximum  extent  and  which  illustrates  to  the  reader  how 
to  structure  a  difficult  program  sequence.   Nine  subroutines 
are  entered  (I,  C,  III,  SI,  S2 ,  S3,  S4 ,  S5,  and  IV)  from  a 
sequence  of  seven  CODE  files  (ITOIV,  S5,  S4,  S3,  S2,  Si, 
ATOD)  on  disk  THESIS9:  using  all  seven  available  $USES  state- 
ments.  Note  also  that  as  long  as  no  more  $USES  statements 
are  required,  an  even  greater  number  of  subroutine  calls 
could  be  present  (e.g.,  S5  could  call  I,  II,  and  III  which 
are  in  the  same  unit,  and  Pascal  "segment,"  as  IV). 
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Table    10:      Contents   of   Demonstration    File   THESI59 : COMPILE. TEXT 

JDSES  DI  IN  THESIS9:ITOIV.  CODE  OVERLAY 
5USES  US5  IN  THESIS9-.S5.CODE  OVERLAY 
SUSES  USU  IN  THE SIS9:S4. CODE  OVERLAY 
SUSSS  US3  IN  THESIS9:S3.C3DE  OVERLAY 
SUSES  US2  IN  THESIS9:S2.CDDE  OVERLAY 
$OSES  aSl  IN  THESIS9:S1.C3DE  OVERLAY 
SOSES    UA    IN    THESIS9: ATOD.30DE    OVERLAY 

PROGRAM    3MPILE 

CALL  I 

CALL  C 

CALL  S1 

STOP 

END 


Table    11:      Contents   of   Deaonstration    File   THESIS9: ITOIV.TEXT 

SUBROUTINE    I 

SRITE(*,'  {  A)  ')      'SUBROUTINE    I    ENTERED' 

RETURN 

END 

SUBROUTINE  II 

WRITE(*,'  (  A)  •>   'SUBROUTINE  II  ENTERED' 

RETURN 

END 

SUBROUTINE  III 

WRITE(*,'  (  A)  •)   'SUBROUTINE  III  ENTERED' 

RETURN 

END 

SUBROUTINE  IV 

WRITE(*,'  (  A)  '>   'SUBROUTINE  IV  ENTERED' 

RETURN 

END 
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Table    12:      Contents   of  Demonstration    File   THESIS9: ATOD.TEXT 

SUSES    UI    IN    rHESIS9:ITOIV.  CODE    OVERLAY 
SUBROUTINE    A 

WRITE(*,»  (A)  •)      'SUBROUTINE    A    ENTERED1 
SRITE(*,»  (  A)  ')      'CALLING    SUBROUTINE    I1 
CALL    I 
RETURN 
END 

SUBROUTINE    B 

WRITE(*,«  (  A)  ')      'SUBROUTINE    3    ENTERED' 
IIITE(*#"  (  A)  ')      'CALLING    SUBROUTINE    II' 
CALL    II 
RETURN 
END 

SUBROUTINE    C 

WRITE(*,»  (  A)  ')      'SUBROUTINE   C    ENTERED' 
WRITE(*,»  (  A)  ')      'CALLING    SUBROUTINE    III* 
CALL    III 
RETURN 
END 

SUBROUTINE    D 

»RITE(*,'  (A)  •)      'SUBROUTINE    D    ENTERED' 
IRITE(*#»  (A)  •)      'CALLING    SUBROUTINE    IV 
CALL    IV 
RETURN 
END 


Table    13:      Contents    cf    Daa onstration    File   THESIS9 :S 1 .TEXT 

SUSES    UI    IN    THESIS9: ITOIV. CODE    OVERLAY 
SUSES    US5    IN    THESIS9  :S5.  CODE    OVERLAY 
SUSES    USU    IN    THESIS9:SU.C3DE    OVERLAY 
SUSES    US3    IN    TBSSIS9 :S3.C0DE    OVERLAY 
SUSES    US2    IN    THESIS9:  S2.C0DE    OVERLAY 

SUBROUTINE    S1 

»RITE(*,«  (A)  ')      'SUBROUTINE    S1    ENTERED' 

»RITE(*,'  (  A)  ')      "CALLING    SUBROUTINE    S2' 

CALL    S2 

RETURN 

END 
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Table   14:      Contents   of   Demonstration.    File   THESIS9:S2.TEXT 

SOSES    01    IN    rHESIS9:ITOIV. CODE    OVERLAY 
SOSES    0S5    IN    THESIS9:S5.CDDE    OVERLAY 
SOSES    0S4    IN    THESIS9:SU.CDDE    OVERLAY 
SOSES    0S3    IN    THESIS9:S3.C3DE    OVERLAY 

SOBROOTINE    S2 

IHITEC*,"  {  A)  ')      'SOBROOTINE    S2    ENTERED1 

WRITE(*,'  (  A)  ♦>      'CALLING    SOBROailNE    S3* 

CALL    S3 

RETURN 

END 


Table    15:      Contents    of   Dam  on strati on    File    THSSIS9:S3.TEXT 

SOSES    01    IN   rHESIS9:ITOIV.CODE    OVERLAY 
SOSES    OS5    IN    THSSIS9:S5.C0DS    OVERLAY 
SOSES    0S4    IN    T3ESIS9:S4.CDDE    OVERLAY 

SOBROOTINE    S3 

WRITER,'  (  A)  •)      'SOBRDOTINE    S3    ENTERED* 

WRITE  (*,'  (  A)  ')      'CALLING    SOBROOTINE    S4» 

CALL    SU 

RETORN 

END 


Table    16:      Consents    of   Dea onstratio*    File    THESIS9:SU.TEXT 

SOSES    01    IN    THES IS9: ITOIV. CODE    OVERLAY 
*OSSS    0S5    IN    TRESIS9 :S5. CDDE    OVERLAY 

S03R00TIME    S4 

WRITE  (*,'  (A)  ')      'SOBROOTINE    S4    ENTERED' 

WRITE(*,'  (  A)  ')      'CALLING    SOBROOTINE    S5« 

CALL    S5 

RETORN 

END 
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The  reader  who  is  developing  a  large  program  for  the 
Apple-II  Plus  is  encouraged  to  try  placing  the  "skeleton" 
of  the  anticipated  program  on  the  machine  before  proceeding 
with  final  coding.  This  will  confirm  that  the  control  flow 
and  overlay  structure  are  possible,  but  will  not  prove  that 
the  final  program  will  fin  into  core  and  leave  room  for  the 
desired  number  of  variables  during  execution. 

Table  17:   contents  of  Daionstration  File  THESIS9:S5.TEXT 

SOSES  01  IN  THESIS9: IT0I7. CODE  OVERLAY 
S0BRO0TINE  S5 

WRITE(*,'  (A)  •)   'S0BR30TINE  S5  ENTERED' 
iRITE(*,'  (  A)  •)   'CALLING  S0BR00TINE  IV 
CALL  IV 
RETOBN 
END 


If  it  becomes  necessary,  during  program  development  or  as  a 
feature  of  the  final  code,  to  know  the  amount  of  free  memory 
available  for  more  instructions  or  variables,  a  Pascal  pro- 
gram must  be  used.   The  reason  for  this  is  that  the  Apple-II 
Plus  is  really  just  a  host  for  the  "P-mcahine"  or  pseudo 
machine  described  in  Reference  [10]  and  this  virtual  machine 
uses  Pascal  accessible  registers  to  keep  track  of  memory 
allocation.   Appendix  B  of  the  above  reference  contains  a 
memory  map  which  shows  that  the  amount  of  free  memory  is 
defined  to  be  the  difference  between  "KP",  which  is  the  top 
of  the  program  stack  (this  translates  to  the  highest  address 
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Chapter  15  of  Reference  [9]  attempts  to  describe  use  of  a 
Pascal  FUNCTION  unit  by  an  Apple  FORTRAN  main  program  and 
dives  an  example  Pascal  program  (complete  with  two  minor, 
unintentional  syntax  errors). 

In  an  attempt  to  make  this  requirement  as  painless  as 
possible,  the  author  has  written  and  tested  a  Pascal  program 
that  will  provide  a  FORTRAN  function  called  MEMREMO  which 
the  reader  may  use  to  find  out  the  amount  of  free  memory 
remaining.   The  two  tables  below  provide  listings  of  the 
Pascal  unit  and  a  FORTRAN  main  program  which  will  illustrate 
the  use  of  the  MEMREM()  function.   In  order  to  use  the  func- 
tion and  test  program  the  reader  must: 

1.  Boot  the  system  with  the  APPLEO:  and  APPLE1: 
diskettes  in  volumes  #4:  and  #5  . 

2.  Enter  the  E(ditor  and  Knsert  the  Pascal  unit  as 
shown  in  the  table  below.   Q(uit  and  U(Pdate  as  usual. 

3.  C(ompile  the  resulting  TEXT  file  using  the  Pascal 
compiler.   The  resulting  CODE  file  need  not  be  L(inked. 

4.  Enter  the  F(iler  and  S(ave  the  TEXT  and  CODE  files 
to  the  desired  FORTRAN  working  diskette. 

5.  Create  the  FORTRAN  program  in  the  usual  manner  and 
C(ompile  using  the  FORTRAN  compiler  after  re-booting 
with  F0RT1:  and  F0RT2:  in  volumes  #4:  and  #5. 

6.  Enter  the  L(inker  and  link  in  the  usual  manner  with 
the  FORTRAN  CODE  file  as  host  and  SYSTEM . LIBRARY  and 
the  Pascal  CODE  files  as  library  files. 
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Table  18:   Contents  of    Pascal  File  #1 1 : MEMREM. TEXT 


(*$P*) 

(♦$S*,L    CONSOLS:    ♦  ) 
UNIT    MEMORY; 
INTERFACE 
FUNCTION    MEMRSM:     INTEGER; 
IMPLEMENTATION 
FUNCTION    MEMREM; 
BEGIN 

MEMREM:=MSMAVAIL; 
END; 
BEGIN 

(♦    THIS    FUNCTION    RETURNS    THE    AMOUST    OF    FREE    MEMORY    ♦  ) 
END. 


Table    19:       Apple    FORTRAN   Prograa    Example    Using   the    MEMREM  () 
Funct  ion 


SUSES    MEMORY    IN    • 1 1 : M EMREM . CDDE 
C 

PROGRAM    HEHTSI 
C 

IMPLICIT     PSAL(A-H,3-Z,     INTEGER  (I-N) 

INTEGER    WORDS, REALS, BYTES 
C 

WORDS    =    MEMREM      () 

3YTES    =    2    ♦    WORDS 

REALS    =    WORDS    /    2 

WRITE  (♦,10)     3YTES,  WORDS,  REALS 

STOP 

10    FORMATf//,'     SPACE    REMAINING    FOR:',/, 
1    4X,1I3,'     BYTES     OF    A    PROGRAM    OR',/, 
1    4X,1I3,»     INTEGER    7ARIABLE3    OR',/, 
1    4X,1l8,f     REAL     VARIA3LES    !!!•) 

END 


C 
C 

C 
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For  comparison  purposes,  the  test  program  given  above 
resulted  in  a  linked  FORTRAN  CODE  file  of  27  (512  byte) 
blocks  of  diskette  storage  (the  bulk  of  which  was  consumed 
by  the  Run  Time  Unit  code)  .   When  executed  the  program 
returned  a  value  of  24512  bytes  of  program  available  out 
of  an  approximate  theoretical  maximum,  prior  to  program 
loading,  of  41,000  bytes.   This  is  enough  room  for  12256 
integer  variables  or  6128  floating  point  variables. 

The  reader  should  also  note  on  the  memory  map  in 
Reference  [10]  that  both  the  program  and  the  variables  in 
the  data  heap  appear  to  be  able  to  write  over  the  two  pages 
of  High-Resolution  Graphics  memory  located  between  8K  and 
24K  in  the  Apple-II  Plus  RAM.   The  author  has  not-  inves- 
tigated the  possible  ramifications  of  this  storage  scheme, 
but  doubts  the  validity  of  the  MEMREM()  function  if  high 
resolution  graphics  are  used  in  the  same  program  with  this 
f unc  _  ion . 


APPENDIX  B 
BLOCK  SOLVER  PROGRAM  LISTING  =  HEWLETT  PACKARD  ENHANCED  BASIC 
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APPENDIX    C 
SAMPLE    EQUATION    SOLUTION    =    HEWLETT    PACKARD    SYSTEM    45 
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APPENDIX  E 
SAMPLE  EQUATION  SOLUTION  =  APPLE-II  PLUS  PERSONAL  COMPUTER 
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